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MULTISENSOR  DETECTION  STUDY 


1.0  INTRODUCTION 

The  intent  of  the  work  which  is  reported  herein  is  to  explore  the 
possibility  that  improvements  can  be  made  in  the  initial  detection  of 
signals  in  noise  using  multiple  sensors,  arbitrarily  located,  by  basing  the 
detection  strategy  upon  a  multivariate  statistical  approach. 

1.1  Background 

This  question  of  detection  is  part  of  the  larger  task  depicted  by 
Figure  1-1,  in  which  the  data  from  m  sensors  are  somehow  to  be  combined 
and  operated  upon  to  produce  "decisions  and  numbers".  The  "decisions"  we 
shall  be  concerned  with  here  consist  of  detection  of  a  source  of  interest 
in  the  medium  being  considered  (e.g.,  underwater)  along  with  whatever 
parameter  estimates  ("numbers"),  are  required  to  carry  out  detection.  We 
restrict  our  attention  to  a  single  source,  whose  waveform  we  denote  by  s(t;e) 
to  indicate  variation  in  time  and  dependence  upon  certain  parameters  0. 

Whether  this  source  of  interest  is  present  or  not,  the  medium  is  such  that 
there  exists  at  each  of  m  sensors  a  noise  waveform  n.ft,^.),  1=1,  2,...,m; 
the  noise  parameters  (n^)  are  in  general  different  in  value  at  each  sensor. 

In  this  work,  the  signal  and  noise  parameters  are  considered  to  be  unknown 
a  priori  and  therefore  must  be  estimated.  Also  both  random  and  deterministic 
signals  are  considered. 

By  "sensor"  we  shall  refer  in  this  work  to  whatever  appropriate  trans¬ 
ducer  and  conditioning  may  be  required  to  acquire  data,  plus  additional 
processing  such  as  sampling  and  analog-to-digital  conversion.  In  some  instances 
a  discrete  Fourier  transform  (DFT)  is  performed  to  obtain  what  is  taken 
to  be  "sensor  data".  The  physical  locations  of  the  sensors  are  assumed  to  be 
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IGUREi-1  'MULTISENSOR  SYSTEM  MODEL 


different,  so  that  observation  in  both  space  and  time  is  performed  by 
collection  of  sensors,  whose  outputs  are  assumed  to  be  available  to  a 
centralized  processor,  not  necessarily  in  "real  time".  Our  primary  concern 
is  with  the  structure  of  this  processor;  therefore,  our  modelling  effort 
begins  with  the  "data"  x^t)  from  the  m  channels,  as  called  out  in  Figure 
1-1. 

1.1.1  Concept  of  Data  Vectors  ,Matrix.. 

At  a  given  sampling  time  tk  the  data  from  the  sensors  may  be 
considered  asan(mxl)  vector  x^.  Over  an  observation  period  if  ji  such 
vector  samples  are  collected,  the  totality  of  the  data  may  be  represented 
by  an  (mxn)  matrix  X  as  illustrated  in  Figure  1-2.  The  columns  of  the 
matrix  are  the  vector  data  samples: 

X  =  {x.j »  Xg*  •  •  '  (1-1) 

The  elements  of  the  data  matrix  are 


xik  ~  x.(tk),i  1,  2,...,m, 

k  =  1,  2,...,n.  (1-2) 

As  discussed  in  Section  1.3,  for  narrowband  signals  and  noise,  the 
{xikJ  can  be  considered  complex  numbers  describing  the  waveform  components 
referenced  to  the  center  frequency.  For  broadband  signals  and  noise,  the 
{xik}  are  more  appropriately  taken  to  be  the  complex  Fourier  coefficients 
of  the  input  waveforms. 

1.1.2  Objectives  of  the  Study 

An  earlier  work  [1]  explored  the  general  subject  of  maximum  likelihood 
(ML)  signal  detection  and  estimation  when  the  data  vector  and  matrix  formulati 
are  used.  It  was  shown  that  in  general,  ML  processors  of  vector  data  are 
analogous  to  their  familiar  scalar  counterparts,  but  are  not  simple 


"building  block"  or  modular  extensions  of  scalar  processors.  For  example, 
the  analogy  to  noise  power  (variance)  on  a  single  sensor  is,  for  multiple 
sensors,  the  determinant  of  the  matrix  of  covariance  between  pairs  of 
sensors.  This  quantity  is  a  nonlinear  combination  of  the  single  channel 
variances  and  also  of  correlations;  thus  something  more  than  just  putting 
together  m  scalar  processors  is  required  to  process  m  sensors  properly. 

In  this  study,  the  theme  of  ML  detection  using  data  from  multiple 
sensors  is  pursued  further,  with  the  following  objectives: 

(a)  Describe  the  theoretical  performance  of  the  ML  detector  of 
multisensor  (multivariate)  signals  in  noise,  using  conventional  figures 
of  merit. 

(b)  Develop  an  computer  program  which  implements  the  ML  detector 
and  which  can  be  used  in  experiments,  also  giving  an  indication  of  the 
complexity  of  the  ML  multivariate  detector. 

1.2  Methodology 

In  accordance  with  its  objectives,  this  study  employs  the  following 
analytical  procedures. 

1.2.1  Theoretical  Base 

In  order  to  determine  the  maximum  likelihood  detection  strategy 
using  the  sensor  data  these  procedures  are  followed: 

(a)  Model  the  joint  probability  density  function  (pdf)  for  the  data 

X  at  the  sensor  outputs  under  the  hypotheses  Hq.-  noise  only,  and  H^:  signal 
plus  noise. 

(b)  Take  the  pdf's  Pj(X;e  ,rj|Hj)  and  pq(X;o|Hq)  to  be  likelihood 
functions  with  e  and  n  representing  signal  and  noise  parameters,  respectively 
Maximize  these  likelihood  functions  by  deriving  maximum  likelihood  (ML) 
estimates  e  and  fi  for  these  parameters,  which  are  NOT  ASSUMED  TO  BE  KNOWN 

a  priori . 
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(c)  Use  the  ML  detection  criterion 


A(X)  =  p^x-.e.n  HjJ/PQCxjnjHQ) 


(1-3) 


to  find  the  statistic  or  function  of  the  data  z( X)  which  tests  the 
hypothesis  Hq.  That  is,  the  above  equation  reduces  to 

z(X)  <  zn  (1-4) 

H0 

where  z^  is  a  threshold. 

(d)  Determine  the  distribution  of  the  test  statistic  z(X)  so  that 
the  probability  integrals 


Q(z0;n|H0)  =  Pr{z>z0|H0}  =  Pr{HQ  rejected  |HQ}  (1-5) 
and 

Q(z0;e,n|H1)  =  Pr{z>z0|Hj}  =  Pr{H0  rejected  |Hj)  (1-6) 

may  be  computed. 

1.2.2  Application 

The  ML  detection  strategy  is  applied  to  the  sensor  detection  problem 
by  carrying  out  the  next  steps: 

(a)  Compute  the  thresholds  zQp  which  satisfy 

Q(z„F;nlH0)  ■  pfa  (i-7) 

for  various  values  of  the  false  alarm  probability  PpA  and  using  those 
thresholds,  compute  the  detection  probability 

Pq  =  Q(zQp;e ,n | H j)  (1-8) 

as  a  function  of  the  actual  values  of  the  parameters  (thus  obtaining 
receiver  operating  characteristics  or  ROC). 


(b)  Relate  the  statistical  detection  parameters  to  the  system 
variables  (number  of  sensors,  number  of  samples)  as  well  as  to  conventional 
input  data  parameters  (signal-to-noise  ratios,  signal  and  noise  inter-sensor 
correlations,  signal  phases,  etc.)  in  order  to  show  the  dependencies  and 

to  permit  comparisons  with  other  systems. 

(c)  Develop  a  computer  program  which  accepts  data  in  the  form 
studied  and  which  performs  the  statistical  procedures  involved  in  the 
detection  strategies;  exercise  the  program  with  simulated  data  to  the 
extent  necessary  to  verify  that  it  works  as  predicted  by  theory. 

1.2.3  Assessment 

The  performance  of  the  multivariate  detection  strategy  is  to  be 
assessed  with  respect  to  the  following  considerations: 

(a)  Minimum  detectable  signal. 

(b)  Complexity  and  storage  requirements  of  implementation. 

(c)  Compatibility  of  configuration  with  other  signal  processing  tasks. 

1.3  •  Models  and  Notation 

A  glossary  of  mathematical  symbols  and  notation  employed  is  given  at  the 
end  of  the  report.  The  concept  of  denoting  the  n  samples  of  waveforms 
at  m  sensors  by  a  data  matrix  X  has  been  introduced  previously.  Although  through 
out  the  report  the  data  is  considered  as  narrowband  (complex)  time  domain  sam¬ 
pled  data,  as  described  next,  the  model  used  can  also  be  extended  to  broadband 
data  as  developed  by  discrete  Fourier  transforms,  as  shown  in  Section  1.3.2. 
Complex  data  are  treated  to  cover  the  subject  thoroughly;  real  (baseband)  data 
is  somewhat  simpler  than  the  complex,  and  can  be  derived  from  it. 

1.3.1  Narrowband  Signals  ana  Noise 

The  notation  which  is  used  is  based  upon  the  assumption  that  the 
vector  )((t)  of  the  m  waveforms,  when  referenced  to  a  given  frequency  and 
phase,  can  be  represented  by  the  narrowband  (Rician)  decomposition 


x.(t)  =  u(t)  cos(u>t  +  4>)  -  j/(t)  sin(ut  +  «f).  (1-9) 

in  which  u(t)  and  v^t)  are  the  in-phase  and  quadrature  components  of  x{t) 
with  respect  to  cos(wt  +  $).  We  may  just  as  well  represent  x(t)  as  the 
(lowpass)  complex  vector  waveform 

x{t)  =  u(t)  +  jv(t),  j  =  v^T,  (1-10) 


(1-11) 


and  the  matrix  of  samples  is 

X=  ||xik||  =  ||uik  +  jv.k||  =  U  +  jV. 

The  columns  of  these  matrices  retain  the  interpretation 

(1-12) 

while  the  rows  are  the  observations  of  the  output  of  a  single  sensor  over 
time: 


(xikW 


~  2L.  =  .£.(  )  * 

0  ^0  Ko 


(x..  )  .  .  =  {x.  (t.),  k  =  1,  2 . n).  (1-13) 

IK  1  lg  »Q  K 

An  example  of  how  this  data  may  be  collected  in  practice  is  diagrammed  in 
Figure  l-3a. 

1.3.2  Multivariate  Gaussian  Model 

'  Now  if  the  waveforms  are  from  stationary,  jointly  Gaussian  random 
processes  with  (mxm)  covariance  matrix  I  =||o.r||  and  mean  vector 
jj  =  ai  +  jl>,  then  the  pdf  for  a  single  vector  data  sample  is 


1 )  =  C ( 2ir )m J e |  D  *  exp|-|{(uk-i)'z  ^(u^-a.) 


+  (yk-b)‘z"1(yk-b)]|. 


(1-14) 


By  noticing  that  the  scalar  quantities  in  the  exponent  may  be  understood 
as  traces  of  matrices  we  can  write 


(u^-a)’!'1^-!)  =  trUu^-a)'!-1^-!)]  (scalar) 
=  tr[E~1(uk-a)(uk-a)'] 


(1-15 
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(a)  NARROWBAND  QUADRATURE  PROCESSING 

Data:  x^k  =  uik  +jvik;  i  =  sensor  index, 

k  =  time  index. 
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(b)  DISCRETE  FOURIER  TRANSFORM  PROCESSING 

Data:  5ik  =  €ikR  +  j£ikI;  1  =  sensor  index. 


complex  Fourier 
coefficients 


FIGURE  1-3.  SCHEMATIC  REPRESENTATION  OF  POSSIBLE  SOURCES  OF  DATA 
AS  MODELLED  (not  necessarily  in  real  time). 
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since  tr(PQ)  =  tr(QP).  Therefore  the  pdf  of  x.  can  alternately  be  written 


P(VH.Z)  =  [(2it)m|E|]"1  etr j-i  (i-ig) 


1  _-l. 


This  form  takes  advantage  of  the  facts  that 

tr[PQ]  =  tr[P-Re{Q}],  P  =  P'(real),  Q  =  Q*; 

trCPQj]  +  tr[PQ2]  =  tr[P(Q1  +  Q2)]; 
and  that  (u^-a )  ( u^-a ) '  +  ( v^-b)  ( v^-b)  ‘  =  Refx^*}. 


(1-17) 

(1-18) 

(1-19) 


For  independent  data  vector  samples,  the  joint  pdf  for  all  the 
observed  data  is  given  by 

P(X)  =  TT  pfx*)  =  [(2n)m| j;j3~n  exp  -ly^C^-aj'E^^-a) 

k=l 


+  (vk-b)*E-A(vk-,b)]j 


II 

=  [(27r)m|z|]"n  etr  -\z~l  .  (2_j 

k=l 

Also,  by  defining  a  matrix  (m  x  n)  whose  columns  all  are  equal  to  the 


mean  vector. 


M  =  (£» 

then  we  may  utilize  the  observation  that 

=  ( X— M)  ( X-M)  ^ 


to  write,  finally. 


(1-21) 


(1-22) 


p(X)  .  t(2.)”|t|r"  etrj-j  r1(y.-M)(X-M)* 


(1-23) 


The  identification  of  the  mean  vector  and  covariance  matrix  with  signal 
and  noise  parameters  will  be  made  in  the  text  in  various  ways.  However 
throughout  the  report  we  shall  assume  that  Gaussian  form  is  suitable  for 


representing  the  data. 


■vlMi 


1.3.3  Extension  to  Broadband  Signals  and  Noise 

In  many  applications  the  narrowband  representation  just  discussed 
is  not  suitable.  However,  the  same  pdf  model  (1-23)  can  be  used  if  the 
interpretation  of  the  data  matrix  X  is  changed  from  time  samples  to  Fourier 
coefficients,  as  illustrated  in  Figure  l-3b. 

Consider  that  the  time  sampling  takes  place  as  previously  described 
and  that  the  n  (real)  samples  from  each  sensor  are  transformed  by  a  discrete 
Fourier  transform  (OFT).  A  new  (complex)  data  matrix  F  =  llc^ll  then 
created.  Its  element  at  the  i:th  sensor  and  r:th  frequency  bin  is  given 
by  [7] 

*ik  » -  H-M-  (1-24> 


Assuming  that  the  time  samples  have  zero  means,  the  covariance 
between  elements  Is 

n  n 


E^ircps(  =  n2 


Ex 


k=r  i= 


cikV  w 


r(k-l)-sU-l) 


If  the  time  sample  vectors  x^  are  independent,  then 

eMps|  •  7  £  Ehkvl  w(r's)tk-1) 

■7  “IP  £“(r'S,(k-l,'^*rs- 

n  K=1 


(1-25) 


(1-26) 


That  is,  if  the  x^  are  independent,  then  so  also  are  the  frequency  sample 
vectors  which  form  the  columns  of  F,  and  it  is  also  evident  from  (1-26) 
that 


11 


i 


L 


c, I 


=  iy\  i 

nZ-J  n 


E^4' 


x^  independent.* 


(1-27) 


For  _x( t)  real,  the  element  obey  the  symmetry  relationships 


’ir 


=  i 


-j  ^  —  1.2 . n-1 


(1-28) 


so  that  in  total  for  each  sensor  the  n  time  samples  produce  n  spectral 
samples  (£  +  1  real,  J  -  1  imaginary).  If  a  subset  of  nf<  J  -  1  frequencies 
is  selected  for  processing  (r=rQ>0  to  r0  +  np  -  1  <  J) ,  then  there  are  np 

complex  vectors  in  the  data  base. 

1.3.4  Distribution  of  OFT  Components 

The  real  and  imaginary  components  ^  and  of  these  vectors  are 

independent,  with  covariance  matrix 


Ej^rR^rR{  Ei^rl^rl\  2n  1 


(1-29) 


Therefore  the  nc  vectors  have  the  joint  pdf 
h  ,  n 


p(F;nF)  *  [(2.)m  exp 


-i 


C*  z  L 


r=l 


(1-30) 


If  the  time  samples  contain  a  deterministic  signal,  then  (1-33) 
is  modified  to  become 


p(F;nF)  =  [(2*)'”  zj]  F  expj-n^^  -  j^)* 


(1-31) 


where  ^  is  used  to  write  the  OFT  component  vectors  of  the  signal. 

It  is  evident,  then  that  solution  of  the  narrowband  detection 
problem  for  the  time-sampling  model  of  the  data  X  will  furnish  also  the 
solution  of  the  broadband  detection  problems  for  the  DFT  model  of  the  data 


♦The  independent  sampling  assumption  is  not  difficult  to  justify  in  the 
case  of  one  sensor  since  it  is  a  matter  of  bandwidth.  For  multiple  sensors, 
one  can  expect  cross  correlations  to  occur  between  different  sensors  at 
different  times  for  directional  sources.  By  ignoring  this,  in  effect  we 
are  requiring  that  either  the  sensor  field  is  small  or  that  the  data  are 
aligned  in  time  (i.e.,  the  array  is  steered  in  the  direction  of  the  source). 
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F  when  it  is  assumed  that  This  assumption  corresponds  to  a  signal 

whose  spectrum  is  "flat"  over  the  bandwidth  spanned  by  the  np  frequency 
bins.  For  this  reason,  only  the  narrowband,  time  domain  case  is  treated 
in  detail. 


m 
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2.0  MULTIVARIATE  MAXIMUM  LIKELIHOOD  DETECTION  STATISTICS 


In  this  chapter  the  likelihood  functions  corresponding  to  hypotheses 
concerning  signal  and  noise  parameters  are  formed  and  used  to  derive  statistics 
or  functions  of  the  data  for  testing  the  hypotheses  (detection).  The  true 
values  of .the  parameters  are  treated  as  unknown  and  are  estimated  so  as  to 
maximize  the  likelihood  functions  under  the  given  hypothesis. 

2.1  Unknown  Deterministic  Signals 

If  the  signal  portion  {s^}  of  the  data  vectors  may  be  considered  to  be 
constant  over  the  interval  in  which  the  samples  are  observed,  then  we  may 
use  the  model 

EO^}  =  jj.  =  JL»  =  !»  k  =  1,  2,...,n.  (2-1) 

This  model  corresponds  to  a  signal  whose  amplitude  and  phase  with  respect  to 
the  reference  frequency  and  phase  are  "slowly  varying"  or  constant  (though 
different  in  value  at  each  sensor),  an  idealistic  extension  of  the  narrowband 
assumption 


Bandwidth  (signal)  «  center  frequency  (2-2) 

The  covariance  matrix  e  of  the  data  in  this  case  corresponds  to  that  of  the 
noise  component. 

If  instead  of  time  samples  the  data  are  Fourier  coeffi cents,  the 
assumption  that  the  mean  vector  is  a  constant  (over  frequency)  corresponds  to 
a  signal  whose  spectrum  is  flat  over  the  n  frequency  bins. 

Detection  of  the  signal  under  these  assumptions  then  is  equivalent  to 
choosing  between  the  hypothesis 


Hq:  v  =  0.  (2-3a) 

and  its  alternative  H ^ ^  f  £.  (2-3b) 

2.1.1  Multivariate  Processing 

The  pdf  of  the  data  given  the  presence  of  a  signal  is 

p(X;£,E|H1)  =  [(2ir )m| E | }~n  etrj4  E'^X-MjtX-M)^.  (2-4) 

Taking  this  function  for  the  likelihood  function  corresponding  to  Hj  and 


maximizing  it  with  respect  to  the  unknown  parameters  and  z  results  in  [1-3] 

i-'sE*  <2-5> 

k=l 

and  „ 

*  =  Tn  5^  (Xfc' £)*  -  ^A-  (2-6) 

k=l 

That  is,  the  ML  estimates  for  the  mean  vector  y^  and  the  covariance  matrix 
z  are  the  sample  mean  vector  y_  and  the  sample  covariance  z.  If  analogous 
to  the  matrix  M  previously  defined,  we  define  the  (mxn)  matrix 

Xq  -  (£,£,  —  »y) »  (2-7) 

then  we  can  write 

i  -  ^  (x-x0)(x-x0)*  •  i  (XX*  -  x0x0*) 

■  ^  (X-M)(X-M)*  -  i  (J  -  u)(C  -  u)*.  (2-8) 

We  are  not  here  concerned  with  the  properties  of  these  estimates,  except 
that  they  are  statistically  independent  as  will  be  discussed.  Substituting 
the  estimates  in  the  pdf  (2-4)  yields  the  Hj  likelihood  function 

Lj(  X)  =  max  p(X;y  ,z  |Hj)  =  [(2Tr)m|£ |  ]'n  e"mn.  (2-9) 


Under  HQ  the  mean  is  assumed  zero  and  the  covariance  ML  estimate  is 


_n 


Vw  =  2n  XX*  =  2n  A0 


(2-10) 


(2-11) 


with  the  resulting  likelihood  function 

L0(X)  =  max  p(x;z|H0)  =  [ (2Tr)m| ZQ | ]_n  e"mn. 

£ 

The  maximum  likelihood  detection  in  this  case  would  be  implemented  by  the 
test  of  the  statistic 


A(X)  = 


IV 


mu 


,A|  "  Txx*  -  x0x0*r  • 


(2-12) 
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the  ratio  of  the  determinants  of  the  estimated  covariance  matrices, 
expression  may  be  reduced  to  another  form  as  follows: 

A(X)  =  |A|-1  | Aq |  =  |A-1|  |A+ n  £  u*|  =  JI  +  nA^  u_  £*1  • 
Middleton  [4]  gives  the  identity 

co 

Jl  +  yGj  =  exp  -  ^ -  tr  Gr  ; 

r=l 

here,  since  tr£(A  1pu*)rj  =  tr  ^-1£(v*A  %)r  V*] 

=  (0*X'1u)r,  (y* A”^y  is  scalar) 


we  have  a(X)  =  exp 


=  exp 


E(-n  v.*A  \ 

r= 1 

|  jl  n  jl  +  n  y*A  ^yj 


y*A  ^w)r/r| 


=  1  +  n  y*A-1y. 


The  ML  detector  performs  the  test 


z(X)  *  ny*A  ^y  5  zn.  (2-1 

H0 

Computational  forms  are  to  be  discussed  in  another  section;  c_r  next 
concern  is  to  determine  the  distribution  of  thi;s  test  statistic  so  that 
probabilities  of  false  alarm  and  detection  may  be  calculated. 

2.1.2  Distribution  of  the  Multivariate  Test  Statistic 

If  the  data  were  real  instead  of  complex,  then  the  test  statistic 
(2-16)  could  be  identified  with  the  "Hotelling  T2  statistic"  [2,  5],  so 


named  because  for  m=l  it  reduces  to  the  square  of  the  familiar  Student's 

t-statistic  (used  for  testing  the  mean  of  a  normal  population  of  unknown 

variance).  It  is  shown  by  Anderson  [2]  that  the  T  statistic  is  distributed 

as  Fisher's  noncentral  F-statistic  with  m  and  n-m  degrees  of  freedom: 

z(X)  distributed  as  —r  Fm  „  _(x) ,  X  real,  (2-17) 

*  n-m  m,n-m  *  *  # 

where  the  noncentrality  parameter  A  is 

A  =  na'  z~la,  (2-18) 

to  use  £ =  Re{£}  to  represent  real  data  (e  is  the  true  noise  covariance 

matrix). 

The  question  now  is,  what  is  the  corresponding  distribution  for 
complex  data?  In  Appendix  A  it  :is  shown  that  in  this  case 

z(X)  distributed  as  JL  F2m,2(n-m)(A>  »  x  complex,  (2-19) 
where  the  noncentrality  parameter  now  is 

A  *  n^'  E"1  a  +  nb'  z~l  b  «  m^E-1  v.  ( 2-^20) 


2.1.3  Combined  Single-Sensor  Processing 

The  stated  objectives  of  this  study  include  determining  what 
improvements  in  detection  performance  may  be  obtained  by  using  multivariate 
processing  of  multiple  sensor  data,  rather  than  some  form  of  modular 
or  "built  up"  approach  based  upon  single-sensor  processing.  The  fact 
that  the  ML  detector  has  been  shown  to  be  indeed  a  multivariate  processor 
is  reason  enough  to  expect  improvements.  However,  it  is  natural  to  wonder 
just  how  much  improvement  can  be  expected,  since  multivariate  processing 
at  least  seems  to  be  a  complicated  procedure. 

In  order  to  provide  a  measure  for  detection  improvements, 'we  now 
consider  ways  of  combining  single-sensor  detection  processing  (performance 
figures  of  merit  are  discussed  in  Chapter  4).  The  test  statistic  for  a 
single-sensor  is  found  by  substituting  m=l  into  the  multivariate  expression, 
(2-12)  to  get  at  the  ith  sensor 


(2-21) 


^  2 

,1,1.'?  I'i  .  “if  +  nl"il 

“i  Xi  -  °0,i/<’l,i - Z~. - 

n 


where 


and 


=  1  +  n|y^|2/a^. , 


■it  \£}hic  “  5it2  s  2n°li 


k=l 


‘1  =  l  X)Xil 

k=l 


(2-22a) 


(2~22b) 


are  the  sample  variance  and  mean.  From  (2-19)  we  know  that 

zi(xi)  *  nly.^/a..  distributed  as  ^  F2,2n-2^V  (2-23) 


and  the  noncentrality  parameter  x;j  is  related  to  the  true  mean  and 
variance  at  sensor  i  by 

Ai  =  n|Ui|2/o2. 


(2-24) 


One  method  for  combining  these  functions  of  the  data  at  single 

sensors  is  to  sum  or  average  them,  as  illustrated  in  Figure  2-1  when 

the  combining  operation  is  addition.  That,  is  we  can  create  a  single 

test  statistic  from  the  m  sensors  by  calculating 

m 

:  *c  ■  <*>  -  s  L  vv 

i=l 

m 

•1+ff  I>il2/“ir  <2-25> 

ial 

Another  method  which  follows  naturally  from  the  multivariate  case  is 
to  combine  the  single  sensor  test  statistic  by  forming  their  product. 

In  effect,  this  is  the  multivariate  solution  for  the  special  case  of 
independence  between  sensor  noises. 


Unfortunately,  the  distribution  for  either  of  these  combinations 
has  not  been  found,  so  their  detection  performances  cannot  be  calculated. 

A  method  which  can  be  evaluated  is  the  following:  perform  individual 
detections  for  each  sensor.  If  the  decision  at  sensor  i  is  "SIGNAL," 
assign  the  value  *1;  if  the  decision  is  "NO  SIGNAL"  at  the  single  sensor 
i,  assign  the  value  y^  =0.  Then  add  up  the  (y^ }  and  decide  for  the 
collection  of  sensors  whether  a  signal  is  present  using  the  rule 


.  SIGNAL 

z~  ?  /  yi  $  k,  0  <  k  <_m.  (2- 

d  "  NO  SIGNAL 

This  type  of  "majority  vote"  technique  may  be  evaluated  fairly  readily. 
The  threshold  (identical)  at  each  sensor  is  chosen  by  requiring  that 


k,  0  <  k  <  m. 


(2-26) 
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Pr(zd  >  k|H0)  =£;(*)Pfr  (l-pF)m-r  =  PFA(given)  (2-27a) 

r=k* 

and  solving  for  Pfr  as  a  function  of  PFA  and  k.  Then  the  required 
individual  sensor  detection  parameters  (all  equal  to  *d)  are  found  by 
requiring  that 

klHl*  =  (r)pn  ^"PD^  =  PD^iven^  (2-27b) 

r=k 

the  solution  for  pQ  can  be  related  to  Xp. 

Without  regard  to  the  distribution  of  the  data  (except  assuming 
sensor  data  are  independent)  we  may  calculate  pF  and  p^  for  given  PFA 
and  PQ.  For  PFA  =  .01  and  PQ  =  .9,  the  results  are 


,i 

Zd 

.9487 

decision  rule 

2  out  of  2  (k  =  2) 

.0589 

.8042 

2  out  of  3  (k  =  2) 

.0420 

.6795 

2  out  of  4  (k  =  2) 

.1056 

.7534 

3  out  of  5  (k  =  3) 

.1504 

.6458 

5  out  of  10  (k  =  5) 

.2498 

.6145 

10  out  of  20  (k  =  10) 

Sensor  data 


VS. 


Sensor  data 


single-sensor  threshold 

test  statistics 
or  decisions 


*  =  some  combining  operation 


FIGURE  2-1.  MULTIVARIATE  VS.  COMBINED  SINGLE-SENSOR  PROCESSING. 
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The  value  of  Pp  determines  the  threshold  and  one  then  computes  the 
X  required  to  achieve  a  single  sensor  detection  probability  of  p^. 

How  this  is  done  is  shown  in  Section  3.1;  each  sensor  decision  is  based 
on  the  F-statistic  distribution  (2-19)  with  m=l. 

2.2  Random  Signals 

If  the  signal  portion  {s^}  of  the  data  vectors  may  be  considered 
to  be  independent  simples  of  an  m-dimensional ,  zero-mean  random  process, 
then  we  may  use  the  model 

UxtfsjL  «  0,  k  =  1,  2 . n.  (2-28) 

The  covariance  matrix  z  of  the  data  in  this  case  corresponds  to 
that  of  the  noise  component  filus  that  of  the  signal  comppnent,  when  present. 

Detection  of  the  signal  under  these  assumptions  then  is  equivalent 
to  choosing  between  the  hypothesis 


H0:  1  ”  Inoise  5  r0  (2-29a) 

and  its  alternative  Hjl  •  Inojse  +  tsigna,  =  ty  (2-29b) 

2.2.1  Multivariate  Processing 

The  pdf  of  the  data  given  the  presence  of  a  signal  is 

ptt-.ElHj)  =  [(21r)m|z1j]'n  etr{-|z"J  XX*},  (2-30) 

where  z,  =  z  +  z  .  Maximizing  with-  respect  to  the  unknown  matrix 

1  noise  signal 

Ij  results  in  the  estimate 

n 

*1  =  In  =  ^  V  *2"31* 

k=l 

and  substituting  the  estimate  in  the  pdf  yields  the  likelihood  function 

l,(X)  .  [(2.)”  |£j|rn  e'™".  (2-32) 

The  Hp  likelihood  function  was  found  previously  (equation  (2-11)). 
Immediately  we  find  that  the  likelihood  ratio  is  identically  one  for  testing 
these  hypotheses  because  both  Eq  and  Zj  are  estimated  by  the  same  quantity. 


» 


I 


> 


I 


I 


» 


Aq/2 n.  Put  simply,  we  cannot  tell  the  difference  between  signal  plus  no^se 
and  noise  only,  under  these  assumptions.  Therefore,  the  hypotheses  need 
to  be  modified  by  additional  assumptions. 

Let  the  noise  at  each  sensor  be  assumed  independent  of  noise  at 
different  sensors,  that  is,  let  us  take  Zq  to  be  diagonal: 


H0:  z  =  dia(oj,  . . .  ,o*)  -  lQ. 
Under  this  hypothesis,  the  pdf  of  the  data  is 


(2-33) 


p{X;  £0IHo)  =  JT  «il  "  e*P  4  542  |xikl2/oi 

t  1-1  J  i=l  k=l 


(2-34) 


Maximizing  this  pdf  with  respect  to  the 


k=l 


and  the  likelihood  function 


2 

noise  powers  (op  yields  the  estimates 

(2-35) 


(2-36) 


The  likelihood  ratio  test  then  becomes 

l1/" .  i£0l  .TnVjj  "i 

J  TTT  lAol  h0 


(2-37) 


that  is,  if  the  product  of  the  diagonals  of  the  sample  covariance  matrix 
Aq  is  greater  than  its  determinant,  then  we  decide  that  there  is  some 
correlation  existing  between  data  from  different  sensors  because  a  signal 
is  present. 

For  example,  suppose  that  we  are  examining  data,  from  two  sensors 
(m=2).  Then  the  test  is 


alla22 


al  l'a  22  — I a  12  ( ' 


1  -  R(x) 


Ho 


(2-38a) 


I 


or 


R(X)  = 


l^l2  "1 
-a  4  '  *  0. 

air 22  Hq 


(2-38b) 


Here  R( X)  is  the  estimated  square  of  the  cross-correlation  coefficient 


between  sensors  1  and  2: 


-  ,  2 


1 23  xikx2kl 


R(X)  a 


(t  i*iki2)gi^i2) 


’12'  * 


(2-38c) 


If  the  data  in  this  case  were  the  Fourier  coefficients  at  a  single 
frequency,  then  R(X)  would  correspond  to  what  is  called  the  "magnitude 
squared  coherence  function"  estimate  at  that  frequency  [8]. 

For  m  sensors  we  can  interpret  the  test  statistic  given  by  (2-37) 
by  further  development: 


1,  (2-39) 

H0 


In  this  relation  we  have  used  (2-14)  and  have  defined  the  matrix  of  sample 
correlation  coefficients  by 

0  p12  p13* ' ' *  plm 

n  A  ^12  0  P23'*‘‘  p2m  .  air 

R  ~  .  o.  =  - ,  i^r. 


plm  p2m  p3m“”  0 


The  trace  of  R  is  zero,  and  we  can  rewrite  (2-39)  as 

H1 

\  trR2  -  4  trR3  +...  S  0. 

C  O  ill 

Ho 


How,  the  traces  of  the  powers  of  R  are 
i  r  i<r 


/a  -  -a 
n  rr 


(2-40) 


(2-41) 


trR3  '2X2}  eir  ^rp  epi>  etc- 


(2-42) 


1  r  p 


Since  | p - r | <1 •  we  can  see  that  an  approximation  to  the  test  (2-41) 
would  be  to  use  the  first  term,  or  the  statistic 


z(X)  *|trR2(X)  =»ir|: 

i<r 

S^tXkxrkl2 


i<r  k 


Kk*  Yz_Jxrkl 


(2-43) 


This  quantity  is  simply  the  sum  of  the  magnitude  squared  of  all  the  measured 
correlation  coefficients  between  sensors. 

2.2.2  Distribution  of  the  Multivariate  Test  Statistic 


In  order  to  set  a  threshold  for  the  test,  we  need  to  determine  its 
distribution.  For  two  sensors  this  task  is  not  too  difficult;  in  Appendix 
D,  the  pdf  for  z(X)  =  |p121  1S  shown  to  be 

p(z)  =  (n-1) (l-pjg)11  (l-z)n_2  2F1(n,n;l;p22z),  (2-44) 


where  p12  is  the  true  correlation  coefficient  and  the  2F1(  )  is  the  Gaussian 

hypergeometric  function.  This  expression  is  in  complete  agreement  with  the 
pdf  for  the  magnitude  squared  coherence  function  estimate  reported  in  [8], 
which  also  gives  the  distribution  function  (see  Figure  2-2) 


P2U0)  s  Pr{z5z0| 

*  2  "]n  n-2  /  .  \  I# 

.  Zfl  _±^_  V  4-^lY 

1_p12z0  k=0  \1“p12z0/ 


The  more  general  case  is  straightforward  but  extremely  complicated 
[9,  10];  so  much  so,  that  the  general  practice  among  statisticians  is  to 
calculate  moments  of  the  test  statistic  (and  these  are  tractable  for  HI  only!) 


(a)  Probability  density  function  {normalised  by  maximum  value) 

for  random  signal  test  statistic,  correlation  varied  (n=32). 


(b)  Cumulative  distribution  function  for  random  signal  test  statistic 
correlation  varied  (n=32). 

FIGURE  2-2.  DISTRIBUTION  OF  RANDOM  SIGNAL  TEST  STATISTIC  FOR  TWO  SENSORS 
(from  [8]) 


Then  the  probability  of  rejecting  Hq  (false  alarm)  is  obtained  by  fitting 
a  pdf  to  the  moments.  In  Appendix  B  it  is  shown  that  the  moments  of 


z(X)  =  |I  +  Rj"1  are 


:/_v\  _TT*  B(n-v-n-k,  k-1) 

r  /  'kU2  B(n+l-k ,  k-1) 


By  noting  that 


B(n-v+l-k,  k-1)  dx  xn-k'v(l-z)k'2 

j  0 
1 

~  f  dx  x"v  pfi(x;n-k+l,  k-1) 

J  0  6 


(2-46) 


(2-47) 


we  are  able  to  say  that  the  inverse  of  z  is  distributed  as  the  product  of 
m-1  independent  beta-distributed  variables.  But  again,  this  is  not  especially 
helpful  for  more  than  a  few  sensors.  So  the  moment  method  is  the  reasonable 
way  to  proceed,  and  is  considered  further  in  the  next  chapter. 
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3.0  THEORETICAL  CALCULATIONS 

A  large  part  of  the  work  consists  of  obtaining  numerical  results 
to  portray  the  theory  developed  and  to  enable  assessment  of  the 
detection  procedures.  In  this  section  the  various  computational 
methods  used  are  explained  and  the  results  tabulated.  Computational 
features  of  any  software  for  implementing  the  detectors  are  considered 
separately  in  another  chapter. 

3.1  Probability  Integral  for  the  F-statistic 

Both  detection  statistics  being  compared  may  be  treated  by 
considering  the  probability 

oo 

Fr{F2n,2(n-m)(l>>4'/dF 

n 

=  Q(n | 2n,2(n-m) ;x)  (3-1) 


This  probability  is  given  by  [6] 

00  k 

Q(n|2n,2n-2m;X)  =  e'X/2  Iy(n-m,  m+k) 


k=0 


where  x  = 


n-m 


n-m  +  mn 

and  Pearson's  incomplete  Beta  function  is  given  by 


Ix(n-m,  m+k)  =  1  -  ^  x(m+k,  n-m) 


■/ 


n-m-1, \m+k-l 
^  B(n-m,~m+k) 


(3-2a) 

(3-2b) 


(3-3) 


A  particularly  simple  computational  form  is 

m+k-1 

Ix(n-m,  m+k)  =  xn"m  ^ 

r^O 


(3-4) 


*» 


using  (n-m)r  =  (n-m)(n-m+l) . . .(n-m  +  r-1) .  (3-5) 

3.1.1  False  Alarm  Threshold 

For  no  signal,  X=0  and  it  is  obvious  from  equation  (3-2)  that 
the  probability  of  false  alarm  is  computed  using 

PFA  =  Q(zol2n‘  2n_2m»0) 


=  I  (n-m,m) , 
*0 


x  _  n-m 
0  "  n-m  +  mZg 


(3-6) 


Since  we  wish  to  calculate  the  threshold  Zg  (or  Xg)  for  given  values 
of  PpA>  it  was  convenient  to  perform  the  calculation  of  (3-6)  using 
(3-4)  on  a  calculator  (hp  34C),  iterating  on  Xg  until  the  desired 
PpA  was  obtained.  (Program  P-1). 

For  the  special  case  of  m=l,  which  applies  to  a  single  sensor, 
(3-6)  reduces  to  the  simple  form  (see  (3-5)) 


Q(z0|2n,  2n-2)  =  xg'1  . 

Table  3-1  gives  zQ  and  xQ  for  PpA 


(3-7) 

.1,  .01,  and  .001  respectively. 


3.1.2  Detection  Probability 

Given  the  threshold  corresponding  to  a  chosen  false  alarm 
probabi li ty  the  probability  of  detection  Pg  is  computed  using  (3-2)  with 
n=  Zg.  That  is,  Pg  =  Q(Zg|2n,  2n-2m,x)  =  Pp(x).  This  was  done  using 
the  FORTRAN  program  P-2.  The  results  are  tabulated  in  Table  3-2. 

The  behavior  of  the  probability  of  detection  can  be  understood 


by  considering  Figures  3-1  and  3-2,  which  are  plots  of  data  selected 
from  the  tables. 


PFA  m 

n 

10 

20 

50 

TOO 

200 

10"1  1 

.774267 

.885867 

.954095 

.977010 

. 988496 

2.6239 

2.4479 

2.3375 

2.3296 

2.3160 

2 

.63164 

.810237 

.922922 

.961279 

.980595 

2.33272 

2.10604 

2.00436 

1.97375 

1.93932 

3 

.50992 

.743495 

.89501 

.947135 

.973477 

2.24255 

1.95499 

1.83779 

1.80471 

1.78913 

4 

.40058 

.68141 

.868716 

.933764 

.966738 

2.24457 

1.87018 

1.73793 

1.70243 

1.68591 

5 

.30097 

.62247 

.84345 

.920878 

.960233 

2.32259 

1.81951 

1.67046 

1.63248 

1.61514 

10 

.35793 

.72559 

.860288 

.929548 

1.79384 

1.51276 

1.46161 

1.44004 

20 

.510666 

.74769 

.872118 

1.43734 

1.34981 

1.31970 

10"2  1 

.599484 

.784760 

.910298 

.954548 

.977124 

6.0129 

5.2112 

4.8285 

4.7140 

4.6589 

2 

.45595 

.6982 

.87207 

.93482 

.96711 

4.77289 

3.89029 

3.52073 

3.41651 

3.33284 

3 

.3437 

.62593 

.83929 

.91779 

.958431 

4.45553 

3.38653 

2.99990 

2.89622 

2.84809 

4 

.25001 

.56125 

.80925 

.90209 

.95041 

4.49976 

3.12695 

2.71069 

2.60488 

2.55670 

5 

.17098 

.50175 

.78095 

.88721 

.94279 

4.84864 

2.97907 

2.52443 

2.41545 

2.36658 

10 

.25395 

.6541 

.81954 

.90792 

2.93778 

2.11527 

1.98177 

1.92695 

20 

.43655 

.69929 

.84515 

1.93603 

1.72009 

1.64900 

_3 

10  J  1 

.464159 

.695193 

.868511 

.932603 

.965883 

10.389912 

8.330548 

7.418375 

7.154455 

7.029047 

2 

.3349 

.60661 

.82662 

.91051 

.95455 

7.943864 

5.836551 

5.033897 

4.815993 

4.713792 

3 

.23885 

. 53459 

.79131 

.89173 

.944863 

7.435699 

4.933357 

4.131727 

3.925773 

3.831945 

4 

.16285 

.47145 

.75936 

.87459 

.936 

7.710930 

4.484463 

3.644332 

3.441430 

3.350427 

5 

.10253 

.41445 

.7296 

.8585 

.92765 

8.753243 

4.238509 

3.335526 

3.131625 

3.041718 

10 

. 18904 

.59904 

.78655 

.890005 

4.289886 

2.677350 

2.442375 

2.348195 

20 

.38378 

.66204 

.82372 

2.408489 

2.041931 

1.926043 

PFA 

=  Q(z0|2m,2n-2m)  =  I 

(n-m,m) , 

Xq  =  (n-m)/(n-m  +  idZq) 

Table  entries: 

xQ  (top) 

*  Z0 

TABLE 

3-1.  FALSE  ALARM  THRESHOLDS  FOR 

F-STATISTIC. 

TABLE  3-2 (a )  PROBABILITY  OF  DETECTION  FOR  F-STATISTIC  (n  =  10,20) 


TABLE  3-2(b)  PROBABILITY  OF  DETECTION  FOR  F-STATISTIC  (n  =  50,100) 


Figure  3-1  shows  PQ  vs  x,  the  noncentrality  parameter,  for  n=50 
and  various  values  of  Pp^  and  m.  This  type  of  plot  is  also  known  as 
"receiver  operating  characteristics"  (ROC).  As  x+0  (-°°  dB),  Pp(x)  -»-Pp^ 
and  this  trend  can  be  observed  on  the  left  side  of  the  figure.  Also 
the  approach  of  Pp(x)  -+  1  as  x  -*•  °°  is  seen  on.  the  right. 

What  is  interesting  in  this  figure  is  that  for  all  PFA* 

Pq(x;  m1)  <  Pq(x;  m2)  for  (3-8) 


This  behavior  is  attributable  to  the  loss  of  degrees  of  freedom  in 
estimating  ^  and  Z;  the  larger  m  given  n,  the  smaller  the  "left  over" 
degrees  of  freedom  2(n-m).  Thus,  for  example,  if  Pp^  =  .1  and  a 
PD  of  .5  is  desired  then  a  larger  value  of  x(6.9  dB)  is  required  for  m=2 
than  for  m=l  (5.6  dB).  At  first  glance,  it  might  seem  that  there  is  a 
disadvantage  in  increasing  m  (more  data  is  bad)!  However,  as  will  be 
discussed  in  Chapter  4,  other  things  being  equal,  the  increased  m 
automatically  insures  a  x  sufficiently  larger  to  produce  an  improved 


V 

Figure  3-2  illustrates  how  for  fixed  m(=3),  the  probability 
of  detection  increases  with  greater  n,  the  number  of  vector  samples 
observed.  This  trait  is  entirely  expected  since  both  degrees  of 
freedom  and  the  amount  of  information  (data)  are  proportional  to  n. 
The  most  informative  aspect  of  the  figure  is  the  demonstration  that 
as  n  increases,  the  curves  converge  to  a  "limit  curve".  This  happens 
because 

PD(x;n,m,z0)  -*■  Q(x2|2m,x),  x2  =  2mzQ  (3-9) 

as  n  -*•  that  is,  the  noncentral  F  distribution  approaches  the 


noncentral  chi-squared  distribution  [  6], 


Therefore,  we  may  write 


Probability  of  Detection 


[•I?  I*  1*1 
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Figure  3-1.  PROBABILITY  OF  DETECTION 
VS.  DETECTION  PARAMETER  FOR  NUMBER 
OF  SAMPLES  FIXED  (n=50)  AND  NUMBER 
OF  SENSORS  VARIED  (m). 

PpA  =  .001,  .01,  .1 
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Figure  3-2.  PROBABILITY  OF  DETECTION 
VS.  DETECTION  PARAMETER  FOR  NUMBER 
OF  SENSORS  FIXED  (m=3)  AND  NUMBER 
OF  SAMPLES  VARIED  (n). 

PFA  =  .001,  .01,  .1 
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3.2  Asymptotic  Distribution  of  Random  Signal  Test  Statistic 

In  Appendix  B  it  is  shown  that  the  moments  of  the  statistic  given 
by  (2-37)  for  testing  the  hypothesis  Hq  (that  the  data  is  uncorrelated 
between  sensors)  can  be  written 


m-1 


e/zv)  =  TT  -UnMjvydd 

1  f 


(3-12) 


Therefore  the  characteristic  function  for  the  new  variable 
y  =  2gn  Jin  z,  0$3<1 

is  given  by 

*y(t)  =  -  E{z32tl!"} 

m-1 


=  JJ  r(n)rf  Bn(l-2jt)-k  +  n(l-e)]  . 

k=l  r(n-k)r[en(l-2jt)  +  n(l-e)] 


We  may  use  the  asymptotic  expansion  [2,  p.  204] 

■i  i 

*nr(x  +  h)  ~  ^  in  2n  +  (x  +  h  -  j)  Jinx  -  x 


L-l 


-E  x'r  wh> +  °<x'L)  <3-i5> 


-Is 


r=l 


(3-13) 


(3-14) 


which  x  and  h  is  bounded  and  the  B  (•)  are  Bernoulli  polynomials 


y  ’ 


[6,  ch.  23].  Applying  this  expansion  to  (3-14)  yields  the  expression 

L 


♦. ,(t)  ~  (l-2it) 


<+\-a/2 


exp 


where  a  =  2  >  k  =  m(m-l) 
k=l 


E 

r=l 


u>r[(l-2jt)‘r  -  1] 


and 


-  (-l)r 


m-1 

(sn)'rE 

k=l 


Br+1[n(l-e)]  -  Br+1[n(l-B)-k] 


(3- 16a) 

(3-16b) 

.  (3-16c) 


Since  (l-2it)“v^  is  the  characteristic  function  of  a  chi-squared  variable 

with  v  degrees  of  freedom,  and  w  «  -y,  the  cumulative  distribution  function 

nr 

of  y  may  be  calculated  using  the  asymptotic  series 

Asymptotic  Order 


Pr|y  $  uj  = 


+  “l(Pa+2  '  Pa} 

+U2+  2  “l)Pa+4““lPa+2"(w2*i“l)Pa 

13  12 

+(u)3+u>1u.2+  6“l) Pa+6’“l(“2+  2  “1 }  Pa+4 

"  wl(w2"iU)l)Pa+2'(a)3'a>lU)2+iwl)Pa 


1 

f1 

-2 


where  P  =  P„ 
v  r 


( 3- 17a ) 
(3-17b) 


This  expression  can  be  simplified  considerably  by  choosing  a  value 
for  the  arbitrary  coefficient  g  in  (3-13)  so  as  to  make  uij=0.  This  value 
turns  out  to  be 


_  n-(m+l)/3 


(3-18) 
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With  this  value  the  asymptotic  series  becomes 

Pr{ysu}  =  Pa  +  »2(Pa+4  -  Pa)  +  »3(Pa+6  -  Pa) 

+[“4<Pa*8  -  Pa>  +  7  “l<Pa+8  '  2Pa*4  +  Pa>] 

+  [“5*Pa+10  ‘  pa*  +  “2“3^pa+10  '  Pa+6  '  Pa+4  +  Paj] 

+  [“6(Pa+12  '  Pa*  +  “2“4(Pa+12  "  Pa+8  "  Pa+4  "  Pa' 

+  7  “3lPa+12  ‘  2Pa+6  +  Pa) 

*  I  “2(PaH2  -  3Pa+8  +  3Pa+4  ’  Pa5]  +  <3'19) 


To  calculate  the  u>r's,  we  need  the  Bernoulli  polynomials;  these  are  listed 
conveniently  in  Table  23.1  of  [6].  After  the  necessary  algebra,  we  have 


“3 

“4 


m-2)(m-l)m(m+l)(2m-l 

162Q[n-(m+l)/3]3 


(m-2)(m-l)m(m+l)(nr-m-7) 

2160[n-(m+l)/3]4 


(3-20) 


Values  of  these  coefficients  for  various  numbers  of  sensors  m  and 
of  samples  n  are  given  in  Table  3-3. 

We  observe  from  these  expressions  that  u)r=0  for  m=2.  This  is 
true  because  [6] 

Br+l(1)  '  W0)  =  lW0>[1-(-1)r] 
and 

B2k+i  h  0,  K  >  0.  (3-21) 

Therefore  for  m=2,  8=(n-l)/n  and  a=2,  with  the  results 

p 

y2  =  2(n-l)tn  z2  is  x2- 


(3-22) 


1 

.004438 

.02400 

.07813 

9.566(-4) 

.004959 

.01543 

.4123 

1.407(-4) 

7. 134(-4) 

.002170 

,05124 

1.079 

3 . 424( -5) 

1 . 724( -4) 

5 . 206 ( -4) 

.01185 

.2307 

8. 446 (-6) 

4 . 237 ( -5) 

1 . 275( -4) 

.002854 

.05356 

B 

1.138{-4) 

1 . 139 ( - 5) 

6 . 426 ( -7 ) 

7 . 712 ( -8) 

9 . 447 ( -9) 

8.960(-4) 

8.415{-5) 

4.592(-6) 

5 . 453( -7 ) 

6.646(-8) 

.003906 

3  -  429 ( -4 ) 

1 . 808 ( - 5 ) 

2. 125( -6) 

2  -  577 ( -7 ) 

10 

.02132 

9.339(-4) 

1 . 039( -4) 

1 . 227 ( -5) 

20 

.04349 

.004299 

4.810(-4) 

B 

- 1 . 969( -6 ) 

-9 . 151 ( -8) 

- 1 . 981 ( -9 ) 

-1 . 172( -10) 

-7.133(-12) 

5 . 760( -5) 

2 . 459 ( -6) 

5.090(-8) 

2 . 971 ( -9) 

1 . 795( -10) 

B 

5.290(-4) 

2.064(-5) 

4.082(-7) 

2.349(-8) 

1 . 410 ( -9) 

10 

.004276 

6 . 604( -5) 

3 . 534( -6) 

2.048(-7) 

20 

.007255 

3.316(-4) 

1.7 88 ( - 5 ) 

Table  3-3  COEFFICIENTS  FOR  ASYMPTOTIC 


EXPANSION 


For  m=3,  for  example,  the  cumulative  probability  distribution 
function  (3-19)  becomes 

Pr{y  £  U}  =  Pr{z2<"-4/3)s  eu} 

*  Pr{4  5  u}  +  2[Pr{*10  5  “}  '  Pr(4  5  “}] 

♦  an^31-3[Pr{42S  Pr{x2  s  u}] 

*  tn-4/3]-4^Pr{x24£u}-  ipr{Xl2su}  *  yy  Pr{x|£u}J 

+  0(n"5).  (3-23) 

2 

The  required  x  probabilities  can  be  found  in  tables  such  as 
Table  26.8  of  [6],  or  calculated  using 

4l-H-£e-U'2  (3-34) 

k=0 

as  shown  by  program  P-4  listed  in  the  back  of  the  report. 

3.2.1  False  Alarm  Probability 

We  now  wish  to  use  the  asymptotic  distribution  (3-19)  to  calculate 
the  false  alarm  probability  for  the  multisensor  problem  with  random  signals. 
The  relation  needed  is 

PFA  ’  Pr{z  >zo} 

■  1  -  Pr{zszo) 

=  1  -  Pr{z2"8  s  22"8} 

=  1  -  Prjy  $  2neJin  zQJ;  (3-25) 

that  is,  the  threshold  zQ  is  found  by  converting  the  value  of  u=uQ  for  which 
Prjy  s  uJ  =  1-Pen  by  the  operation 


(3-26) 


For  large  numbers  of  samples  it  is  evident  from  (3-16b)  (3-18),  and  (3-20) 
that  to  an  excellent  approximation, 

20  =  exp  (xl-PFA,m(m-l)^n  '  *  .  (3‘27) 

Also,  since  z(x)  =  |I  +R|”*»expj|-  trR2|  as  shown  in  Section  2.2.1,  we  may 

interpret  (3-27)  as  a  statement  about  the  average  values  of  the  measured 

1  2 

intersensor  correlations.  That  is  z > zQ  is  equivalent  to  ^  trR  >  £nzQ,  or 

m^2 ^average  >  ^Xl-PpA»  m(m-l)^n"  ^3‘28^ 

Both  Zq  and  this  ! p I average  deflned  by  (3-28)  are  given  for  n=100,  200 
and  several  values  of  m  in  Table  3-4.  For  example,  for  m=5  sensors  and 
n=l 00,  the  ratio  of  the  determinants  of  the  estimated  covariance  matrix 
under  Hq  and  Hj  respectively  has  a  1%  chance  of  exceeding  the  value  1.4672 
when  Hq  is  true;  this  value  would  be  exceeded  by  an  average  inter-sensor 
correlation  coefficient  of  about  0.2  if  were  true. 


*  a 


4.60517 


1.0476 


'average 


.2157 


n  =  200 


1.0234 


'average 


.1521 


10.6446 


1.1139 


1.0550 


.1336 


18.5494 


1.2076 


.1773 


1.0980 


28.4120 


1.3363 


.1703 


1.1543 


.1198 


10  90 


107.565 


3.0544 


1 . 7296 


20  380 


415.728 


87.3735 


.1534 


8.6195 


.1065 


9.21034 


1.0975 


1.0474 


.2151 


16.8119 


1.1858 


.2383 


1.0883 


26.2170 


1.3055 


.2108 


1.1413 


.1484 


37.5662 


1.4672 


.1958 


1.2089 


.1377 


10  90 


124.116 


3.6270 


1.8817 


20  380 


447.1 


122.4285 


.1591 


10.1409 


.1104 


10  90 

20  380 


13.816 

22.458 

32.909 

45.315 

137.208 

470.9 


1.1498 

1.2556 

1.3975 

1.5877 

4.1550 

158.1337 


1.0719 

1.1197 

1.1805 

1.2572 

2.0114 

11.4719 


.2635 

.1941 

.1663 

.1513 


.1133 


Calculation  of  Multivariate  Test  Statistic 
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The  statistics  z(X)  derived  in  Chapter  two  were  derived,  for 
analytical  purposes,  from  forms  originally  expressed  as  ratios  of  the 
determinants  of  sample  covariance  matrices.  Computationally  it  is  simpler 
to  calculate  these  determinants  than  to  perform  the  matrix  inversion  called 
for  by  the  analytical  form.  In  this  section,  the  particular  method  used 
for  computing  determinants  is  documented,  and  it  is  shown  how  to  convert 
thresholds  for  the  F-statistic  to  those  for  the  ratio  of  determinants. 

3.3.1  Computing  the  Determinants  of  the  Sample  Covariance  Matrices 

A  well  known  iterative  method  for  computing  the  determinant  of  a 
symmetric  matrix  is  the  "forward  Doolittle,  left  side"  scheme  [15].  That 
basic  approach  is  adapted  here  for  the  complex,  Hermitian  sample  covariance 
matrices  of  the  data.  It  is  to  successively  reduce  the  dimension  of  the 
determinant  by  subtracting  a  multiple  of  the  first  row  from  the  other  rows 
so  as  to  make  the  first  column  all  zeros  except  for  the  "top"  or  "pivot" 
element.  For  example,  for  a  3»3  determinant,  the  algorithm  works  as 
f ol 1 ows : 


all  a12  a13 
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(3-28) 


An  mxm  determinant  is  found  by  an  algorithm  which  may  be  written 
as  follows: 


Let  D  =  | A | 

Initially,  D  =  DQ  =  1 

For  i  =  1  to  m.  -1 


3.3.2  Computer  Implementation 

The  application  of  this  algorithm  to  the  complex-valued  sample 
covariance  matrices  involved  in  the  detection  problems  under  consideration 
was  tested  by  means  of  the  simple  simulation  listed  as  program  P-5. 

In  this  program,  a  Gaussian  random  number  generator  (ostensibly  producing 
independent,  zero-mean,  unit  variance  numbers)  was  used  to  develop  five¬ 
dimensional  complex  data  vectors,  x^. 


The  mean  vector  and  Hq  sample  covariance  matrix  were  "built"  or 
accumulated  iteratively  by  the  relations  (n  =  number  of  samples  *  10) 


[k] 


=  yi  ^-U  X 


ik 


/n 


l,2,...m(=5) 


(3-29) 


*  ^  rki 

(Eo  £  E  °0M 


=  °0j£  +  Xik‘  V2n) 


i  ,sl  =  1,2,.  ..m. 


(3-30) 


A  diagram  of  the  simulation  is  shown  in  Figure  3-3.  The  estimated 
covariance  matrix  under  the  hypothesis  is  formed  by  subtracting  one 
half  the  outer  product  of  the  mean  vector  from  the  Hq  hypothesis  covariance 
matrix.  Then  the  determinants  of  these  matrices  are  computed  by  the 
algorithm  given  above,  and  the  test  statistic  is  taken  to  be  (the  real  part 
of)  their  ratio. 

Figure  3-4  shows  the  results  of  the  calculations.  The  ratio  for 
ten  samples  equals  1.813  for  this  particular  set  of  random  vectors.  In 
the  next  section  it  is  shown  how  to  calculate  threshold  against  which 
this  number  is  to  be  tested. 

3.3.3  Obtaining  False  Alarm  Threshold  for  Ratio  of  Determinants 

The  statistic  z(X)  for  which  the  false  alarm  detection  probabilities 
were  computed  in  Section  3.1  related  to  the  ratio  of  covariance  determinants 


—  =  R(  X)  =  1  +  nil*  E_1 


=  1  +  z(X) 


1  +  111  C 

1  n-m  r2m,2(n-m) * 


(3-31) 


FIGURE  3-3.  DIAGRAM  OF  MULTIVARIATE  TEST  STATISTIC  CALCULATION. 
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FIGURE  3-4  SAMPLE  OUTPUT  OF  PROGRAM  P- 
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Therefore  R 

>  Rq  implies 

F  *  T  <"b  -  »  •  Fo- 

(3-32a) 

or 

R  >  1  +  7T5  F0  '  R0' 

(3-32b) 

Now  from  Table  3-1 

we  find  for  n=10  and  m=5  that 

the  following  values  may 

be  determined: 

!ea 

Fg(from  table) 

h 

.1 

2.32 

2.16 

.01 

4.85 

3.43 

.001 

8.75 

5.38 

Thus  the  test  value 

,  being  less  than  Rq  in  each 

case  results  in  the  (correct) 

acceptance  of  Hq  at 

the  10%,  1%,  and  0.1%  levels 

• 

4.0  APPLICATION  OF  THEORETICAL  RESULTS 


Having  learned  how  to  calculate  the  probability  of  detection 
PQ  for  the  multivariate  detection  approach  we  are  now  in  a  position 
to  apply  these  results  to  the  detection  problem.  In  this  chapter, 
the  objective  is  to  find  reasonable  values  for  parameters  which  describe 
the  detector  performance  and  to  show  how  the  multivariate  detector  may 
be  implemented. 

4.1  Detection  Performance  Parameter,  Deterministic  Signals 

It  is  convenient  to  extract  a  single  parameter  from  the  tables 
and  curves  of  Chapter  3,  one  which  represents  the  performance  of 
detection  (likelihood  test)  in  achieving  a  useful  Pp  while  rejecting 
false  alarms.  In  receiver  operating  characteristics  (ROC)  for  single 
sensor  detectors,  the  noncentrality  parameter  A  which  we  have  employed 
becomes  signal-to-noise  ratio  (SNR).  It  is  common  to  describe  the 
detector's  performance  by  specifying  the  SNR  for  which  given  Pp  and  PFA 
are  obtained.  For  example,  we  might  choose  the  pair 

(PD.  PFA)  =  (-9.  IQ'2)  '  (4-D 

and  define  the  "minimum  detectable  signal"  (MDS)  as 

MDS  =  SNR  such  that  Pp(SNR)  =  .9  for  PpA  =  10~2.  (4-2) 

and  abbreviate  this  statement  by  the  notation 

MDS( .9)  h  SNR( . 9, .01) (or  SNP,(.9)  with  PpA  =  .01  understood) 

(4-3) 

Since  the  multivariate  parameter  a  is  analogous  to  SNR,  we  shall 
extract  values  of  a  for  which  Pp  =  .5  and  .9  from  the  ROC  data  already 
computed,  then,  in  the  subsequent  sections,  interpret  this  parameter 
in  terms  of  minimum  detectable  signal  by  employing  certain  assumptions 
about  signals  and  noise  at  the  sensors. 


As  an  example,  from  Figure  3-1  we  see  that  for  P 

m  =  2,  a  A  of  6.9  dB  is  required  to  produce  Pp  = 

to  give  Pp  =  .9.  This  we  write  shorthand  as 

A ( .  5)  =  6.8  dB  P...  =  .1 

for  hA 

a( .9)  =  11.3  dB  n=50,  m=2. 


^  =  .1,  n  =  50,  and 
.5  and  a  A  of  11.3  dB 


(4-4) 


Table  4-1  gives  A(.5)  and  x(.9)  values  extracted  from  the 
data  generated  in  Chapter  3.  The  same  information  is  plotted  in 
Figures  4-1,  4-2,  and  4-3.  From  these  curves  we  see  that  A ( . 5)  and 
A(  .9)  decrease  to  an  asymptotic  value  with  increasing  n,  and  increase 
with  m  (the  number  of  sensors)  for  a  given  n.  This  is  in  keeping 
with  the  discussion  made  of  Figures  3-1  and  3-2  in  the  last  chapter. 
4.1.1  Structure  of  the  Multivariate  Detection  Parameter 

From  Chapter  2  (or  Appendix  A)  we  recall  that  the  multivariate 
noncentrality  parameter  is 

A  =  n  u*  £  *  _u  (4-5) 

where  £  is  the  mean  of  the  complex  data  vectors  {x^}  and  l  is  the 
covariance  matrix  of  their  independent  real  and  imaginary  parts.  What 
signal  and  noise  models  correspond  to  these  quantities? 

In  Chapter  1,  the  complex  data  matrix  elements  were  identified 


with  the  in-phase  and  quadrature  components  of  the  received  sensor  wave¬ 
forms  with  respect  to  a  center  frequency: 
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gure  4-1.  REQUIRED  DETECTION  PARAMETER  VALUES  vs.  NUMBER 
OF  SAMPLES,  Pcfl  =  .1 


FIGURE  4-2.  REQUIRED  DETECTION  PARAMETER  VALUES  vs.  NUMBER 
OF  SAMPLES,  Pctt  =  .01 


Figure  4-3.  REQUIRED  DETECTION  PARAMETER  VALUES  vs.  NUMBER 
OF  SAMPLES  Pc.  =  .001 
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or  Xl(tk)  *  ujk  cos„ctk  -  vjk  sfn„ctk 

-1 

cos  [to  t.  +  tan  1  — —  ].  (4-6b) 

C  K  U.k 

If  we  are  considering  only  the  possibility  that  a  deterministic 
signal  is  present  in  the  background  noise,  then  the  covariance  matrix 
Z  *1 |o-r| |  is  the  covariance  matrix  of  the  zero-mean  noise  at  the 
sensors: 

°ir  =  E{ninr}*  (4_7) 

The  mean  vector  p.  =  <L  +  j.b  components  then  are  due  entirely  to  the 
signal,  when  present,  and  we  assume  that  the  signal  components  (in- 
phase  and  quadrature  amplitudes  or,  alternatively,  envelope  and  phase) 
are  constant  (or  nearly  so)  while  the  n  samples  are  observed.  Let 
the  signals  at  the  sensor  outputs  have  envelope  values  { S - }  and  phases 
{ e - } ;  then  the  elements  of  the  mean  vector  can  be  identified  as 

Uj  =  a..  +  jb^  =  S.  cose.  +  jS^  sine^  (4-8) 

Now,  using  oir  to  denote  elements  of  the  inverse  of  i,  we  can 
express  the  multivariate  noncentrality  parameter  by 

m  m 

A  =  n  p*  z"1  v  =  n  (ai  -  j^)  o1r(ar  +  jbp) 

i=l  r=l 

■  “"^Vr  +  fjV  +  j(afbr  -  a^)] 

i,r 

=  n  o1r  S.SrCcos(e.-  er)-  j  sin(e.-er)].  (4-9) 
i,r 

Since  z  is  a  symmetric  matrix,  so  also  is  z'*,  and  the  summation  over 
the  imaginary  part  is  zero,  leaving 


» 


m  m 


c  C0S(B.j  ®r) 

i=l  r=l 


(4-10) 


The  subscript  "M"  stands  for  "multivariate". 

For  the  special  case  in  which  the  noise  is  independent  from 
sensor  to  sensor*  (the  covariance  matrix  is  diagonal),  (4-10)  reduces 


m 

v»Z>"  s< 


Xs  sf 

=  2n  y  - 2“  =  2mn  x  (average  SNR).  (4-11) 

i=l  2oi 

With  this  result  as  motivation,  we  choose  to  define  minimum 
detectable  signal  for  the  multivariate  detector  by 

MDW  ■  25i  W  •-  5S  xCP0^,«1.  (4-12) 


4.1.2  Performance  Predictions 

Using  the  definitions  for  minimum  detectable  signal  that  have 
been  given,  we  may  predict  the  performance  of  detection  based  on  the 
test  statistics.  For  the  multivariate  detector,  we  use  the  numbers  in 
Table  4-1  to  calculate 

MDSM(PD)(dB)  =  \(PD;  m,  n)(dB)  -  10  log1Q(2nm).  (4-13) 

Similarly,  for  the  detector  based  on  a  single-sensor  detector  output 
the  MDS  is 

MDS^PgHdB)  =  A(P0;  1,  n)(dB)  -  10  log10(2n).  (4-14) 

For  example,  for  a  false  alarm  probability  of  .001  and  a  detection 
probability  of  .5,  for  m=10  sensors  and  n=50  samples,  we  have 


MDSm(.5)  =  15.4  dB  -  10  log1()(l000)  =  -14.6  dB  (4-15) 
and  for  m  =  1 , 

MDSJ(.5)  =  11.4  dB  -  10  log1Q(100)  =  -8.6  dB.  (4-16) 


Thus  in  this  example,  the  multivariate  detector  has  a  processing  gain 
of  6  dB  over  a  single-sensor. 

We  may  define  a  "multivariate  processing  gain"  by 


mpgs!!W 

"  MKm(Pd)  ■ 

From  (4-12)  we  may  then  write 


(4-17) 


MPG  «m 


X(PD;  1,  n) 

X( Ppi  m,  n)  ' 


(4-18) 


Using  the  numbers  for  Aj(Pg)  and  x^(Pg)  given  in  Table  4-1  the  MPG  is 
calculated  and  plotted  in  Figures  4-4,  4-5,  and  4-6  for  the  various 
Pp^  values. 
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Figure  4-4.  MULTIVARIATE  PROCESSING  GAIN  vs.  NUMBER  OF  SAMPLES, 
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MULTIVARIATE  PROCESSING  GAIN  vs.  NUMBER  OF  SAMPLES 


It  is  evident  from  these  figures  that  detection  based  upon  the 
multivariate  approach  realizes  appreciable  gains  over  using  a  single 


sensor.  Gains  from  2  to  9  d  B  are  shown.  The  curves  rise  to  an 
asymptotic  or  limiting  value  as  n,  the  number  of  observations,  increases 
and  this  value  is  directly  proportional  to  m,  the  number  of  sensors. 

On  each  figure  we  see  that  the  MPG  is  higher  for  Pp  =  .9  than  for 
Pp  =  .5,  with  the  amount  of  the  increase  a  function  of  Pp^:  the  more 
stringent  (smaller)  the  PpA  required,  the  smaller  the  increase.  By 
comparing  the  figures  we  notice  that  MPG  values  tend  to  increase  as 
Pp^  decreases.  This  behavior  is  very  interesting  because  it  implies 
that  the  multvariate  detector  achieves  gain  in  two  dimensions  (Pp,  PpA) i 
the  more  difficult  the  detection  criterion  (either  Pp  or  PpA)>  the 


better  the  multivariate  approach  is  compared  to  single-sensor  processing. 
That  is,  this  gain  can  be  obtained  for  independent  sensor  noise  terms. 

If  the  off-diagonal  components  of  the  inverse  covariance  matrix  e”*  are 
nonzero,  actually  the  is  more  of  less  than  2mn  •  SNR  (see  (4-10)),  and 
the  performance  of  the  multivariate  approach  can  be  slightly  better  or 
worse  than  has  been  shown,  depending  on  the  amount  of  correlation  between 
sensor  waveforms,  and  the  signal  phases  (see  Appendix  C). 

Comparison  with  Majority  Decision  Scheme.  Analysis  of  the  "majority 


decision"  method  described  in  Section  2.1.3  results  in  the  following 
graphical  comparison: 


m,  number  of  sensors 

Figure  4-7  GAIN  OF  MULTIVARIATE  DETECTOR  OVER  MAJORITY  SCHEME 
These  results  show  the  combination  of  single-sensor  detector 
decisions  (which  assumes  independent  sensors)  doing  better  than  the 
multivariate  for  small  numbers  of  samples  n  or  for  large  number  of 
sensors  m.  However,  the  comparison  is  not  ideal  since  the  multivariate 
processor  is  "working  harder"  than  necessary:  even  though  the  sensor 
noises  are  independent  (in  this  comparison),  it  is  using  up  degrees  of 
freedom  (see  Section  3.1.2)  estimating  correlations.  It  seems  reasonable 
to  expect  that  (a)  a  multivariate  processor  derived  for  the  special  case 
of  independent  sensors  will  uniformly  do  better  that  the  combination 
scheme,  and  (b)  when  inter-sensor  correlations  are  introduced,  the 
performance  of  the  majority  decison  scheme  will  deteriorate,  while  the 
multivariate  will  not.  However,  the  distributions  in  either  of  these 


cases  have  not  been  found,  so  these  expectations  remain  conjectures. 

4.2  Computer  Program  for  Performing  Multivariate  Detection  of 
Deterministic  Signals 

In  order  to  assess  the  implementation  costs  that  may  be  associated 
with  the  multivariate  detection  processing,  a  computer  program  in  FORTRAN 
was  written  and  exercised.  Limited  simulations  using  the  programs  provided: 
verification  of  theoretical  detection  performance  predictions.  j 

„ 


4.2.1  Program  Requirements  and  Desirable  Features 

The  basic  requirement  the  program  must  satisfy  has  already  been 
indicated  in  Section  3.3;  the  program  must  accept  n  vectors  {x^}  which 
are  m  *  1,  accumulate  mean  vector  and  sample  covariances,  calculate  determinants, 
and  form  test  statistics  which  are  ratios  of  determinants.  The  basic 
operation  is  specified  by  the  equation 

l2n£ol 

2n  fQ  -  n 

Although  the  detection  problem  has  been  formulated  id  this  study  as  a  "batch" 
decision  rather  than  a  sequential  or  an  iterative  one,  certain  inherent 
relationships  in  the  functions  to  be  computed  allow  accumulation  of  intermediate 
results.  Therefore,  only  a  minimum  of  input  data  has  to  be  stored..  Specifically, 
the  sample  mean  vector  can  be  accumulated  using 


-  S  zn.  (4-19) 

p.  p*  0 


z(X|  ,x_2  *  • 


■V  = 


where 


-  [n-1] 
p 


+  V" 


-[in] 

R 


+  xi+1/n, 


(4-20) 


Similiarly,  the  Hq  covariance  can  be  accumulated  using 


where 


+  y  y  *• 

-i +1-1+1  ' 


0. 


Certain  economies  are  also  possible  due  to  the  fact  that  the  covariance 
matrices  are  positive  definite  and  Henmitian. 

In  practical  detection  systems  it  is  desirable  to  have  some  control 
over  "integration  time"  or,  in  this  context,  the  number  of  samples.  For 
weak,  stationary  signals  one  should  like  to  make  the  number  of  samples  high; 
for  strong,  intermittent  signals  a  relatively  short  observation  time  would 
be  selectable.  Implementation  of  this  second  desirable  feature  is  a  simple 
matter  of  storing  a  table  of  thresholds. 

How  to  accommodate  the  variable  integration  time  feature  is  suggested 
by  Figure  4-8.  A  100-sample  data  base,  for  example,  can  be  processed  as 
a  whole  (maximum  integration  time)  or  as  ten  10-sample  data  bases  in  succession 
(minimum  integration  time).  In  fact,  both  (plus  in-between  cases)  can  be 
done  in  parallel,  and  the  results  of  each  displayed  simultaneously  in  some 
fashion. 

Assuming  simultaneous  10/20/50/100-sample  detections,  the  basic 
requirements  for  the  multivariate  detector  computer  implementation  may  be 
summarized  as  in  Figure  4-9. 

4.2.2  Program  Structure  and  Size 

The  basic  requirements  for  computer  implementation  of  the  multivariate 
detector  sunmarized  in  Figure  4-9  are  met  by  the  FORTRAN  program  listed  as 
program  P-6  in  the  back  of  this  report.  A  flow  diagram  is  presented  in  Figure 
4-10.  The  symbology  used  is  as  follows: 


samples 


STORAGE :  (not  including  that  required  for  arithmetic,  subroutines,  overhead) 
DATA:  one  m  x  1  input  buffer  for  (complex)  data  vector 

□ 

MEAN  VECTORS:  four  m  *  1  (complex)  arrays 

□  □  □  P 

COVARIANCES:  eight  m  *  m  (complex)  arrays 


ARITHMETIC:  complex  conjugation,  complex  and  scalar  arithmetic 

SUBROUTINES: 


Compl ex 

Vector 

Zeroing 

Complex 

Vector 

Addition 

Complex 

Vector 

Outer 

Product 

Compl ex 

Complex 

Complex 

Matrix 

Matrix 

Matrix 

Zeroing 

Addition 

Determinant 

TABLES:  False  alarm  thresholds  vs.  PpA  vs.  number  of  samples 

FIGURE  4-9.  SUMMARY  OF  BASIC  REQUIREMENTS  FOR  MULTIVARIATE  DETECTOR 
COMPUTER  IMPLEMENTATION. 


START 


STOP 


FIGURE  4-10  FLOW  OF  MULTIVARIATE  DETECTOR  PROGRAM 


Program  Variable  Name 
M(=5) 

.  NSAMP 


.Mathematical  Equivalent 

m,  number  of  sensors 

n,  number  of  samples 


CX(M) 

(complex) 

x^,  data  vector 

CVX(M) 

(complex) 

p,  sample  mean  ( 10 * s  base) 

CVXX(M) 

(compl ex) 

p,  sample  mean  ( 50 ' s  base) 

CVL(M) 

(complex) 

p,  sample  mean  (20’s  base) 

CVC(M) 

(compl ex) 

p,  sample  mean  (100's  base) 

CAX(M,M) 

(complex) 

2n£g,  Hq  covariance  (10's  base) 

CAXX(M,M) 

(complex) 

2n£g,  Hg  covariance  (20 * s  base) 

CAL(M,M) 

(complex) 

2n£0,  Hq  covariance  (50 ’s  base) 

CAC(M,M) 

(compl ex) 

2n£0.  Hq  covariance  (100's  base) 

CBX(M,M) 

(complex) 

2n£j,  Hj  covariance  (10's  base) 

CBXX(M,M) 

(complex) 

2n£j,  Hj  covariance  (20’s  base) 

CBL(M,M) 

(complex) 

2n£j,  Hj  covariance  ( 50 ' s  base) 

CBC(M.M) 

(complex) 

2n£j,  Hj  covariance  (100's  base) 

RATX,  RATXX, 

RATL,  RATC 

z,  computed  test  statistics 

THRX ,  THRXX, 

THRL,  THRC 

Zg,  stored  thresholds 

All  vector  and  matrix  operations  (except  for  construction  of  CVX  and  CAX) 
are  performed  using  subroutines.  The  main  programs  and  subroutines,  compiled 
as  listed  on  an  IBM  3033  FORTRAN  system  were  sized  for  m=5  as  follows: 


Program  or  Subprogram 

Source  Lines 

Size 

MAIN  (incl.  data  generation) 

85 

4442 

CLRV  (zero  vector) 

7 

380 

CLRM  (zero  matrix) 

8 

472 

ADDV  (add  vectors) 

7 

588 

ADDM  (add  matrices) 

9 

734 

PRODV  (vector  product) 

8 

718 

DETER  (determinant) 

15 

910 

complex  arithmetic. 

random  no.  generator 

500 

I/O,  FORTRAN,  other  system  overhead 

8744 

21208 

29952 

The  "size"  is  understood  to  be  in  units  of  8-bit  bytes,  so  that  the 
program  itself  requires  about  8k  on  an  8-bit  machine,  or  4k  on  a  16-bit 
machine. 

As  shown  in  the  flow  diagram,  the  program  operation  is  straightforward: 
10-sample  mean  vectors  CVX  and  Hg  covariance  matrices  CAX  are  constructed 
from  the  input  data  (one  input  vector  at  a  time).  After  incrementing  the 
other  (20- ,  50- ,  and  100-sample)  mean  and  covariance  accumulations,  the 
covariance  CBX  is  formed,  its  determinant  computed,  and  the  test  statistic  z  is 
computed  as  the  ratio  RATX. 

In  lieu  of  decision-making  the  program  prints  out  z  and  the  threshold 
Zg  corresponding  to  a  one  percent  false  alarm  probability. 

FORTRAN  conditional  (IF)  statements  are  used  to  "count"  the  number  of 
10-sample  blocks  of  data  which  have  been  "observed,"  and  appropriately 
enable  the  20-,  50-,  and  100-sample  detection  computations  at  the  proper 
"time".  After  each  detection  operation  the  accumulated  means  and  covariances 


just  used  are  reset. 


4.2.3  Comparison  of  Computed  Multivariate  Results  with  Theor 


To  exercise  the  multivariate  detection  computer  program  and  to 
collect  simulated  results  for  verification  of  the  theory,  a  data  base  of  100 
vectors  with  the  following  characteristics  was  generated  and  processed. 
m=5,  z  ®  I,  ^  = . SJ[  +  j(). 


The  signal  components  were  chosen  to  be  equal  in  phase  (0°)  and  amplitude 
at  each  sensor.  On  separate  runs  the  amplitude  was  made  0,  .2,  .3,  .5,  and 
.9  to  correspond  to  SNR's  of  -17,  -13.5,  -9,  and  -3.9  dB  respectively. 

The  results  of  these  simulations  and  their  theoretical  interpretation 
may  be  displayed  as  shown  in  Figure  4-11.  Detections  for  a  particular 
block  of  data  are  indicated  by  shading  the  portion  of  a  diagram  (similar  to 
Figure  4-9)  which  corresponds  to  that  block  in  simulated  time.  Just  to  the 
right  of  these  diagrams,  a  tally  of  detections  is  provided  (1/5  =  one 
detection  in  five  trials,  etc).  The  other  experimental  values  presented  are 
the  arithmetical  average  values  of  the  test  statistics  compute,  normalized 
by  the  1%  false  alarm  threshold. 

Thus,  for  example,  on  the  run  for  which  a  SNR  of  -13.5  dB  was 
simulated,  the  20's  detection  yielded  a  detection  at  the  third  data  block  only, 
and  the  average  value  of  the  five  statistics  computed  was  83.95%  of  the  false 
alarm  threshold. 

Along  with  the  experimental  data  in  Figure  4-li  theoretical  predictions 
are  given  for  probability  of  detection  and  mean  value  of  the  test  statistic 
relative  to  the  threshold.  The  Pq  values  were  determined  by  setting  x=  2nm  SNR 
and  using  Table  3-2.  The  mean  value  of  the  test  statistic  was  computed  using 


[14] 


v2(vi  +  a) 

"  Vj(v2  -  2) 


(4-22) 


l 


to  write 


E)2(  =  1 


+ 


m 

n-m 


2(n-m)(2m  +  2nm  SNR) 
2m(2n  -  2m  -  2) 


=  1  +  d  +  nSNR)-  (4-23) 

The  Zq's  were  obtained  from  Table  3-1  as  described  in  Section  3.3.3, 
resulting  in 


n 

10 

20 

50  ! 

100 

zQ  (m=5,  PpA  =  1%) 

5.8486 

1.9930 

1.2805  | 

1.1271 

The  experimental  values  are  very  close  to  theoretical  predictions. 
For  example,  for  SNR  =  -9  dB,  the  50-  and  100-sample  detection  calculations 
yielded  "perfect  scores";  theory  predicted  that  ditection  was  better  than 
99%  probable.  The  20-sample  detection  calculations  scored  3/5  or  60% 
compared  to  the  theoretical  65%  probability,  and  the  10%  score  of  the 
10-sample  detector  was  consistent  with  a  predicted  of  12%. 

Another  way  of  viewing  the  results  is  given  in  Figure  4-12, 
in  which  the  theoretical  and  experimental  average  detection  statistics 
are  compared  as  functions  of  SNR.  The  agreement  is  very  good  for  this 
parameter  also,  considering  the  relatively  small  sample  sizes. 
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Figure  4-12  AVERAGE  (NORMALIZED)  DETECTION  STATISTIC  vs.  SNR, 


4.3  Modification  to  Computer  Program  for  Detecting  Random  Signals 
The  computer  program  for  detecting  deterministic  signals  may  be 
easily  modified  in  order  to  calculate  test  statistics  for  random  signals. 
All  that  is  required  is  the  following  subroutine: 


SUBROUTINE  RHODET  (CMAT,  M,  Cl) 

IMPLICIT  COMPLEX  (C) 

C2  =  (1.,  0.) 

D01  I  =  1,  M 

1  C2  =  C2*CMAT (1,1) 

Rl  =  REAL  (C2/C1) 

WRITE  (6,2)  Rl 

2  FORMAT  (IX,  'CORRELATION  DET  RATIO  F9.4) 

RETURN 

END 

In  response  to  the  statement,  for  example  CALL  RHODET  (CAX,  M,  CNX) 
the  subroutine  calculates  the  ratio  of  the  product  of  the  diagonal  elements 
of  the  covariance  matrix  CAX  to  its  previously  computed  determinant,  CNX, 
and  prints  it  out. 


5.0  ASSESSMENTS  AND  RECOMMENDATIONS 


The  study  is  concluded  in  this  present  phase  by  reviewing  what 
has  been  done,  assessing  the  results,  and  recommending  further  work. 

5.1  .Review 

5.1.1  Summary  of  Work  Performed 

The  mathematics  of  complex  multivariate  Gaussian  statistics  have 
been  applied  to  a  model  of  a  multi  sensor  detection  problem.  Testing 
hypotheses  concerning  the  presence  or  absence  of  sources  in  the  sensor 
medium  have  been  shown  to  be  equivalent  to  computing  statistics  or 
functions  under  the  Gaussian  noise  model,  corresponding  to  ratios  of 
determinants  of  sample  covariance  matrices. 

For  the  case  of  deterministic  signals  whose  magnitude  and  phase 
remain  constant  during  the  period  in  which  data  are  collected,  the 
multivariate  test  statistic  has  been  shown  to  be  in  the  family  of  F- 
statistics.  Using  the  properties  of  F-statistics  distributions, 
theoretical  predictions  have  been  made  of  SNR's  required  to  detect  the 
signal  while  rejecting  false  alarms  at  a  given  level.  The  results  confirm 
the  expectation  that  multivariate  processing  of  sensor  data  requires  a 
significantly  smaller  SNR  than  does  single-sensor  processing.  Also,  a 
non-ideal  comparison  indicates  that  multivariate  processing  is  better  than 
a  majority  decision  method  for  combining  individual  sensor  decisions  when 
large  numbers  of  samples  are  used. 

A  computer  program  implementing  the  multivariate  detection  processing 
has  been  developed  and  tested,  yielding  simulation  results  compatible  with 
theory. 

For  the  case  of  random  signals,  the  multivariate  test  statistic  has 
been  shown  to  have  a  distribution  which  approaches  the  chi -squared  family 
asymptotically.  Computation  of  the  random  signal  test  statistic  may  be 
accomplished  by  a  very  slight  modification  to  the  program  for  the  detection 
o f  deterministic  signals. 


5.1.2  Summary  of  Multivariate  Detection  Performance 

For  deterministic  (narrowband  or  broadband)  signals  the  performance 
of  the  multivariate  detector  may  be  summarized  as  follows: 


number 

of 

sensors 


Minimum  Detectable  Signal  in  dB  (Pg 
z  =  diagonal) 
number  of  samples 
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10 

20 

50 

100 

200 

1 

0.6 

-3.1 

-7.4 

-10.5 

-13.6 

2 

-1.0 

1 

cn 

• 

o 

-9.5 

-12.7 

-15.8 

3 

-1.6 

-6.0 

-10.8 

o 

• 

H 

-17.1 

4 

-1.6 

-6.7 

-11.5 

-14.7 

-17.9 

5 

-1.3 

-7.1 

-12.1 

-15.4 

-18.5 

10 

X 

-7.4 

-13.8 

-17.2 

-20.5 

20 

X 

X 

-14.7 

-18.8 

-22.3 

Multivariate  Processing  Gain  in  dB  over 
Single-Sensor  Processing  {Pn  =  .9,  Pcft=.01,  z  =  diagonal) 


2.2 

2.2 

3.5 

3.5 

4.2 

4.3 

4.9 

4.9 

6.7 

6.9 

8.3 

8.7 

5.2  Assessments 

The  multivariate  detection  approach  is  now  assessed  with  regard 
to  detection  performance,  implementation,  and  compatibility  issues. 

5.2.1  Detection  Performance 

Clearly  the  multivariate  detection  approach,  achieving  a  SNR 
processing  gain  of  around  2nm,  is  capable  of  performing  very  sensitive 
detections.  It  is  the  ML  solution  for  the  unknown  signal  case. 

It  is  instructive  to  note  that  a  beamformer  utilizing  m  sensors 
and  n  samples,  and  processing  the  sum  of  m  unknown  signals  in  unknown 
noise  by  the  ML  method  shown  for  single  sensors,  gets  better  results  than 
the  multivariate  technique  when  steered  at  the  target  (equal  signal 
phases).  The  gain  of  the  beamformer  over  the  multivariate  is,  for  the 
parameters  used  in  Section  5.1.2.,  given  below: 

Beamformer  Gain  over  Multivariate  in  dB 


n  = 

10 

20 

50 

100 

200 

m  =  2 

i.4 

1.2 

0.9 

0.8 

0.8 

3 

2.6 

1.9 

1.4 

1.3 

1.3 

4 

3.8 

2.4 

1.9 

1.8 

1.7 

5 

5.1 

3.0 

2.3 

2.1 

2.1 

10 

X 

5.7 

3.6 

3.3 

3.1 

20  . 

X 

X 

5.7 

4.7 

4.3 

Thus  it  is  seen  that  the  multivariate  detector  is  not  far  behind  the 
beamformer  in  performance  when  their  number  of  samples  is  large.  Moreover, 
it  must  be  noted  that  the  multivariate  detector  can  operate  where  good 


beamforming  is  not  possible,  such  as  sensor  placements  which  are  not 
well  known  or  subject  to  random  variation.  This  is  because  in  effect 
the  multivariate  processor  estimates  the  signal  phases  attach  sensor. 

Whether  the  multivariate  processing  gain  is  sufficiently  great 
to  be  worth  implementing  is  primarily  a  function  of  the  effort  required 
to  collect  the  data.  On  the  one  hand,  combined  single-sensor  processing 
or  some  "majority  vote"  of  single-sensor  decisions  seems  a  simple 
procedure;  one  does  not  have  to  push  m  channels  through  the  same  I/O  or 
FFT  device  in  roughly  the  same  time  frame  or  buy  a  parallel  processor. 

On  the  other  hand,  if  the  detection  scenario  tolerates  a  non-real  time 
solution,  such  as  time-sharing  FFT  and  I/O  devices,  the  only  significant 
additional  burden  the  multivariate  approach  imposes  is  the  need  for 
storing  (buffering)  the  data  segments  from  each  sensor  (rows  in  the  X 
matrix)  until  the  entire  data  matrix  is  collected.  In  this  connection 
it  should  be  noted  again  that  the  time  alignment  of  sensor  data  is  not 
critical  to  the  performance  of  the  multivariate  detection  approach,  since 
for  independent  sensor  noises,  the  sensor  signal  relative  phases  do  not 
affect  it. 

The  comments  above  apply  to  the  detection  of  deterministic  signals. 
Evaluation  of  the  detection  of  random  signals  remains  to  be  done. 

5.2.2  Implementation 

As  demonstrated  by  the  computer  program  for  implementing  the  multi¬ 
variate  detector,  the  programming  complexity  and  storage  requirements  are 
quite  modest.  This  assessment,  however,  is  based  upon  the  availability  of 
input  data  vectors  (columns  in  the  X  matrix)  in  sequence.  If  the  data 
collection  situation  differs  from  this  assumption,  additional  storage  for 
buffering  will  be  required.  Still,  the  advent  of  microprocessor-based, 
parallel  processing  at  low  cost  would  make  an  m-channel  data  bus  going  into 
a  multivariate  processor  (another  microprocessor)  a  Strong  possibility. 


5.2.3  Compatibility 

Since  some  detection  scenarios  there  is  a  certain  amount  of 
competition  for  computing  resources  among  the  various  signal  processing 
tasks  (detection,  localization,  tracking,  etc.),  it  is  appropriate 
to  ..assess  the  compatibility  of  the  multivariate  detection  approach  with 
these  other  tasks.  While  it  has  not  been  shown  in  this  work,  in  the 
previous  study  [1]  it  was  demonstrated  that  the  statistical  quantities 
developed  by  the  multivariate  processor  (means,  variances,  and  covariances) 
are  also  sufficient  statistics  for  estimation  of  certain  parameters. 

Thus  simultaneous  detection  and  estimation  can  be  carried  out  using 
some  of  the  same  data  processing. 

For  example,  the  diagonal  elements  of  the  sample  covariance  matrix 
are  estimates  of  noise  power  at  each  sensor  and  the  off-diagonal  elements 
are  inter-sensor  correlations.  The  mean  vector  contains  estimates  of 
the  signals  arriving  at  the  sensors.  Also,  the  single-channel  test 
statistic  is  an  estimate  of  SNR  at  one  sensor,  while  the  multivariate 
test  statistic  is  a  generalized  SNR.  Therefore,  it  appears  that  the 
multivariate  processor  is  very  compatible  with  signal  processing  tasks 
requiring  the  estimation  of  these  quantities. 

5.3  Recommendations 

Further  work  recommended  may  be  classified  into  two  categories: 
continuations  and  extensions. 

5.3.1  Continuations 

Recommended  further  work  in  a  continuation  of  the  present  work: 

(a)  Theory:  Find  receiver  operating  characteristics  for  random 


signals.  Deal  with  DFT  components  explicitly  for  both  random  and 
deterministic  signals.  Obtain  ROC  for  special  cases  of  the  ML  detector 
such  as  independent  sensor  noises. 


(b)  Comparisons:  Compare  multivariate  random  signal  detection 
performance  against  that  of  combined  two-sensor  correlation  detectors. 

(c)  Algorithms:  Refine  the  program  for  computing  test  statistics 
so  that  detection  results  for  the  various  sample  sizes  (10/20/50/100) 
are  available  every  ten  samples  by  using  a  "sliding  window"  concept. 

Also,  seek  ways  to  iterate  the  computation  of  determinants  so  that 
accumulation  and  storage  of  covariance  matrices  is  not  required. 

(d)  Validation:  Test  the  multivariate  detection  approach  against 
simulated  and  actual  multi  sensor  data. 

5.3.2  Extensions 

Recommended  further  work  extending  the  present  work: 

(a)  Develop  joint  detection  and  estimation  procedures  which 
utilizethe  same  block  of  multisensor  data.  For  example,  joint  random 
signal  (correlation)  detection  and  localization  via  time  delay  estimation 
are  especially  compatible.  Also,  a  generalized  multivariate  detector 
could  be  used  to  "educate"  a  beamformer  and  thereby  obtain  bearing 
estimates. 

(b)  Study  in  detail  the  impact  on  software  and  hardware  requirements 
that  implementation  of  multivariate  processing  would  have,  for  a  specific 
multi  sensor  scenario  such  as  airborne  ASW  surveillance  using  sonobuoys. 
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P-1  Incomplete  Beta  Function  (False  Alarm  Probability) 

P-2  F-Statistic:  Probability  of  Detection 

P-3  Asymptotic  Distribution  Coefficients 

P-4  Chi-Square  Cumulative  Distribution 

P-5  Tryout  of  Matrix  and  Vector  Operations 
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PROGRAM  P-1.  INCOMPLETE  BETA  FUNCTION  (F-STATISTIC  FALSE  ALARM 
PROBABILITY). 

Ix(b,a)  =  xb(l  +  Cj  +  c2  +  cfll) 
where  c  ,  =  (b  +  r) (1  -  x)c  /(r  +  1),  cn  =  1.  (a  >  1) 


34C  keystroke  codes: 


LBL  A 

CHS 

STO  I 

RCL  6 

X 

5WT 

STO  4 

RCL  4 

ISG 

GTO  0 

x+-»y 

RCL  0 

RCL  1 

GTO  2 

LBL  1 

STO  0 

1 

STO  6 

GTO  1 

rclT 

RTN 

- 

LBL  0 

LBL  2 

RCL  3 

LBL  B 

EEX 

X 

RCTT 

RCL  1 

STO  3 

3 

ST0+5 

INT 

ytX 

1 

A 

1 

RCL  4 

X 

STO  5 

i 

STO+6 

R+ 

RTN 

- 

+ 

CLx 

A 

Routine  A  stores  values  of  a  and  b. 

To  compute  I,  enter  x,  run  routine  B. 
Example:  Ix(8,2)  *  -10000  for  x  =  .63164. 


PROGRAM  P-2.  PROBABILITY  OF  DETECTION  FOR  F-STATISTIC. 

//D2NFXACT  JOB  ( 0591 , 0000,  ,  ,  ,  ,  , Y, O ) , MILLER, CLASS=E 

/*FETCH 

/^NOTIFY 

/*N03ETUP 

//3TEP1  EXEC  FORTGCLG 
//FORT. SYS IN  DD  * 

DIMENSION  AM ( 20 > ,  PD(20), NN( 10) , MM( 10) , XX (7,  5) 
RE  AD  (5,  1 )  (NN(  I),  1  =  1,5),  (MM(I),  1=1,7) 

1  FORMAT (1215) 

READ (5,  2)  (AIKI),  1  =  1,  11) 

2  FORMAT (11 F5 . 1 ) 

READ (5,  3)  XX 

3  .  FORMAT ( 7F 11.6) 

WRITE <6,  3)  XX 
DO  32  N A=  1,2 
NO=NN ( N A ) 

WRITE <6,  4)  NO 

4  FORMAT! /IX,  'N=',  I3/2X, 'M=  L=0  -6DB  -3DB 

CDB  SDB  10DB  12DB  14DB  16DB" ) 

DO  30  MC=  1,3 
MA=MM(MC) 

MB=NO-MA 
X— XX (MC,  NA) 

(continued  on  next  page) 


PROGRAM  P-2  (continued) 


>o  co 

xO  O'  ^ 

o  o  *-» 

ONN 

in  n  co 


o  o 


o 

CM 


CO  GO 
CO  O'  0 

cn  in  cm  in 

MOOCN 
10  CM  vO  CM 
CO  fv  CO  CN 


O  CO  00 

r*  •  N  N  in  I''  CO 

xo  ex  0  cm 

^  O  01  CO  o  o 

O  CM  «*  CM  '© 

UO  CO  <1  CO  O'  Ox 


I- 

< 

n 


.  x©  0 

CM  CO  -•  -H  x£,  CO 

H  b'l  M  N  N  N 

O  «-•  TO  CO  >0 

(0  O  00  xo  CO  xo 
•  vj  0  Ox  Ox 

o  •  •  •  •  • 


_J 

CM 

IL 

H 

• 

\ 

H 

0 

in  m  N 

* 

CM  O-  rt  CO  N 

tH 

H 

Ox  «a-  o  *«  «* 

1 

ll 

H 

ex  co  in  ex  o 

0 

l-« 

• 

o  «t  ox  «a-  rv 

+ 

O 

^  h 

U0  N  TO  Ox  O' 

CCi 

'*>■  H 

• 

Ox  £ 

Z 

01  >-t 

o 

• 

~  O 

s  ^ 

o 

CO 

o 

0  1- 

1-  1- 

<  P  • 

o  . 

in  <-* 

i-  <r 

<  0 

£  £L  N 

CM  CO 

o 

s 

0  0 

n  n 

<  •«*  Ll 

N  CM  O'  U0 

1- 

n  ^ 

0  -I 

H 

_l 

1  *>  H 

o 

<*  00  CM  N  O' 

o 

H  n 

u. 

O 

lL  "** 

w  ^  *-4 

o 

'O  CM  O'  CM  U0 

o 

n  ^ 

<->  X  * 

#: 

\  UO 

Cl  £  - 

o  • 

O  CM  —*  O 

H 

0  £ 

O  <  -> 

CM 

•  j 

x  <r 

•“*  o 

0  rt  CM  'O  CO 

<■*»  H 

< 

r(  .  £  X 

0  >- 

CM  Ui 

UJ  —  i-i 

xO  0  O'  O'-  o 

■  *. 

^  W 

10  0  1 

£  + 

\ 

)4!  H  H  k 

o  -< 

H  # 

Mill  x  • 

<  . 

-1  £  co  X 

#  o 

•  II 

.  .-i  * 

£  >  H  rt  CM  UJ 

#  £ 

£  I- 

£  0  UJ  -  ti!  UJ 

in  ■ 

O  H  , 

v-  **  • 

+  X  II  «  U  j 

x  0 

<  -J 

~t  ifi—i'G**’  —i 

Q  0 

III 

_  w 

P  £  x©  O 

2<n»:  +  z 

*  CO 

V  •  N 

0  II  Z  -  5-  Z  Z 

p  i 

.  UO  II 

•  »  <C  H  ■ 

•  II  £  CM  CM  w  CM  II  ^ 

H  H 

||  ^  »-i  UJ  <  >-i  >-i 

N  MOO  >0 

X  CC  H  rt  H  II  p  II 
ll 


rt  rt  X  P 

II  II  5  ~ 

.'-.i  Ar  ■ 


m  c  i 
.1!. 


£  II  ~  <  I-  <  O 

«_  g  ?  ^  ij.  £  c«  e  ii  Cm  Cm  £  u.  o  «‘M 
Z  >•  0  0  •-«  P  U 


II  Z  II  £  Z 
CM  O  CM  _»  11 


O  n 
II  ~ 

rl  Q. 


£  -  Z 

_i  ij  p 


>-  O  >  CO  Z  O  M  O  (0  0.  O 


H  £  >-  I-  CL 
m£CZ?OQ 
CC  0  0  U  !-  Z 
3  U.  O  O  CO  UJ 


o  in  ^  o  cm 

^  CM  CO  CO  CO 


'0  *0  O'-  •*  0- 

CM  CO  O  O 

<*  UO  «f  N  CO 

N  CD  in  N  OX 

N  0  O'  O-  O'- 


in  x©  n 


CO  O' 


» 


PROGRAM  P-3.  ASYMPTOTIC  DISTRIBUTION  COEFFICIENTS. 


//D2AXC0EF  JOB  ( 0591 , 0000, , , , , ,  Y,  O) ,  MILLER,  CLASS=E 
/*FETCH 
'  /^NOTIFY 
/*N03ETUP 

//STEP1  EXEC  FORTGCLG 
//FORT. SYSIN  DD  * 

DIMENSION  XM(6),  XN(6),  D2 (6,6),  D3 ( 6,  6>,D4<6,  6) 

READ (5,  1)  <XM(I),  1  =  1,5),  (XN(I>,  1  =  1,5) 

1  FORMAT < 10F7. 2) 

DO  2  NM=1, 5 
DO  3  NN=1, 5 

IF( XM<  NM) .  GE. XN(NN ) )  GOTO  3 
X~XN(NN)—( XM(NM) +1 . )/3. 

D2(NM,  NN)=(XM<NM>-2. ) * < XM(NM) -1 . > *XM (NM) XM(NM) +1 . >/72./X**2 

D3 < NM,  NN)=D2(NM, NN>*<2. *XM(NM>-1 . )/22.5/X 

D4  <  NM,  NN ) =D2 <  NM, NN ) *  <  XM  <  NM  > **2-XM ( NM  > -7. > /30. /X**2 

3  CONTINUE 

2  CONTINUE 

4  FORMAT (6E12. 4) 

WRITE (6,  4)  D2 
WRITE (6,  4)  D3 
WRITE <6,  4)  D4 
STOP 

END 

/* 

//GO. SYSIN  DD  * 

3.  .  4.  5.  10.  20.  10.  20.  50.  100.  200. 


I 


I 


PROGRAM  P-4  CHI-SQUARE  CUMULATIVE  DISTRIBUTION. 

2 

Pr(x  £X  }  =  P(x|\>)  =  1  -  Q(x|v),  where  v  is  the  number  of  degrees  of  freedom 

and  x/2  V/2_1  n 

Q(x| v)  =  e'x/2  l  (x/2)n/n! 

n=0 

hp  34C  keystroke  codes: 


LBL  A 

2 

i 

1SG 

ex 

CLR  £ 

4 

RCL  7 

GTO  1 

STO  9 

EEX 

STO  7 

x-<-+y 

RCL  1 

X 

3 

1 

Jt 

RCL  7 

RCL  9 

J. 

STO  8 

RCL  8 

2 

X 

STO  I 

LBL  1 

X 

4 

x-*-+y 

1+ 

STO  8 

CHS 

Example:  Q(447.1|380)  =  .0100;  Q(4.60517 j 2)  =  .10000. 

PROGRAM  P-5  TRYOUT  OF  MATRIX  AND  VECTOR  OPERATIONS. 


/ /D2TRY0UT  JOB  ( 0591 , OOOO, ,  ,  ,  ,  , Y,  0> , MILLER,  CLASS=E 

/*FETCH 

/^NOTIFY 

/*NOSETUP 

//STEFl  EXEC  FORTGCLG 
//FORT. SYSIN  DD  * 

IMPLICIT  COMPLEX (C) 

DIMENSION  CX<5), CC(5,  5),  CMU(5>,  CC0(5,  5) 

-i  CALL  ST ARTR  ( 65539 ) 

DO  1  K=  1,5 

CMU <  K ) =CMPLX  <  O. ,  0.  ) 

DO  2  Kl=l , 5 

CC(K,  K1)=CMPLX(0.,0. ) 

CCO  <K,  K1 )=CMPLX  <0. ,  O.  > 

2  CONTINUE 

1  CONTINUE 

C  GENERATE  RANDOM  VECTORS,  MEAN  VECTOR,  AND  XX*  MATRIX 
DO  8  N=l, 10 
DO  3  J=1 i 5 
A=RANDN  < 1 . ) 

B=RANDN(1.  ) 

CX  <  J  >  =CMPLX  <  A,  B) 

CMU ( J ) =CMU ( J ) +CX  <  J  > / 1 0 . 

DO  4  J 1 = 1 , J 

CC ( J,  J1 ) =CX ( J  > »CONJG ( CX ( J1 ) ) /20. +CC( J,  J1 ) 

CC  <  J1 ,  J ) =CONJG ( CC ( J,  J1 ) ) 

4  CONTINUE 

3  CONTINUE 

©  CONTINUE 


(Continued  on  next  page) 


PROGRAM  P-5  (continued) 


WRITE ( 6,  12) 

12  FORMAT  ( IX,  '  MEAN  VECTOR  ') 

WRITE (6,  5)  (CMU(I),  1=1,5) 

5  FORMAT (IX,  10F9.3) 

WRITE (6,  6) 

6  FORMAT(//,  IX,  '  MATRIX  XX«  ',/) 

DO  7  1=1.5 

WRITE(  6,  5)  <CC(I,  J),  J=l,5) 

7  CONTINUE 

C  GENERATE  ( X-XO) ( X-XO) *  MATRIX 
DO  13  J=l,5 
DO  14  1=J,  5 

COO ( J,  I ) =CC ( J,  I )-CMU( J)*CONJG(CMU( I ) ) /2. 

CCO ( I,  J) =CON JG ( CCO <  J,  I ) ) 

14  CONTINUE 

13  CONTINUE 
WRITE ( 6,  15) 

15  FORMAT  (//,  IX,  '  MATRIX  <  X  —  X  O  >  <  X  —  X  O  >  ** 
DO  16  1=1,5 

WRITE(6,  5)  (CCOtl,  J),  J=l,5) 

16  CONTINUE 

C  COMPUTE  DETERMINANTS,  RATIO 
CDET=CMPLX ( 1 .  ,  0. ) 

CDETO=CDET 
DO  9  1=1,4 
CPIV=CC(I,1) 

CPIV0=CC0( I, I ) 

CDET=CDET*CPIV 
CDET0=CDET0*CP I VO 
*  UMIN=I+1 
DO  10  J=JMIN, 5 
DO  11  J1=JMIN, 5 

CC(U,  J1)=CC( J, J1 )-CC(I, J1 )*CC< J,  I )/CPIV 
CCO ( J,  J1 ) =CCO ( J,  J 1 ) -CCO  < I , J 1 ) *CCO ( J, I ) /CPI VO 
11  CONTINUE 

lO  CONTINUE 

9  CONTINUE 

CDET=CDET*CC (5,5) 

CDETO=CDETO*CCO (5,5) 

CRAT =CDET /CDETO 
WRITE(6,  17) 

17  FORMAT ( /,  IX,  '  DETS  AND  RATIO  ',/) 
WRITE(6, 5)  CDET, CDETO, CRAT 

STOP 

END 

/* 

//LKED.SYSLIB  DD 

//  DD  DSN=C591 . SUBLIB,  DISP=SHR 

/* 


PROGRAM  P-6  FULL-SCALE  IMPLEMENTATION  OF  MULTIVARIATE  PROCESSING. 


//D2BIG0NE  JOB  (0591, 0000,  ,  ,  ,  ,  ,  Y,  0) ,  MILLER,  CLASS=E 

/»FETCH 

/^NOTIFY 

/«NOSETUP 

//STEPl  EXEC  FORTGCLG 
//FORT. SYS IN  DD  * 

IMPLICIT  COMPLEX <C) 

DIMENSION  CX  ( 5 ) ,  CVX ( 5 ) ,  CAX  ( 5,  5>,CBX<5,  5),  CAXX(5,  5>,  CBXX( 
1 ,  CAL ( 5,  5),  CBL(5,  5),  CACC5,  5),  CBC(5,  5) ,  CVXX  (5) ,  CVL<5) ,  CVC< 
CALL  STARTR( 65539) 

M=5 

N5AMF-10 

ENSAMP=FLOAT  ( NS  AMP  > 

THRX=5. 8486 
THRXX=1 . 9930 
THRL=1 . 2805 
THRC=1. 1271 
WXX=2.*ENSAMP 
WL=5. *ENSAMP 
WC= 1 0 . SENSAMP 
CALL  CLRV(CVXX, M) 

CALL  CLRV(CVL,M> 

CALL  CLRV(CVC,  M) 

CALL  CLRM ( CAX  X , M, M ) 

CALL  CLRM(CBXX, M,  M) 

CALL  CLRM ( CAL, M, M) 

CALL  CLRM ( CBL, M, M ) 

CALL  CLRM ( CAC, M, M ) 

,  CALL  CLRM ( CBC, M,  M ) 

DO  1  1  =  1,5 
DO  1  J=l,2 
CALL  CLRM ( CAX,  M,  M) 

CALL  CLRMtCBX,  M,  M) 

CALL  CLRV(CVX,M> 

DO  2  K- 1 , NSAMP 
DO  2  K1  =  1,M 
A=RANDN< 1. ) +. 3 
B=RANDN ( 1 . ) 

CX ( K1 ) =CMPLX ( A, B) 

CVX ( K1 ) =CVX ( K 1 ) +CX ( K1 ) /ENSAMP 
DO  2  K2= 1 , K 1 

CAX ( K1 , K2 ) =CX ( K1 ) *CON JG ( CX ( K2 ) ) +CAX ( K1 ,  K2 ) 

2  CAX ( K2,  K1 ) =CONJG <  CAX ( K1 , K2 ) > 

CALL  ADDMCCAXX,  CAX, M, 1 . , 1 . ) 

CALL  ADDV (CVXX,  CVX, M,  1.,  .5) 

CALL  ADDM(CAL, CAX, M, 1 . , 1 . ) 

CALL  ADDV (CVL, CVX, M, 1., .2) 

CALL  ABDM(CAC, CAX, M, 1. , 1. ) 

CALL  ADDV (C VC, CVX,  M, 1 . ,  . 1 ) 


(continued  on  next  page) 


cn  cn 


ouo  .  nuuu  in  o  o  o  <t  co  <0^000 


TENS  DETECT 


CALL  PRODV(CVX,  li,  CVX,  M,  CBX,  ENSAMP) 

CALL  ADDM (CBX,  CAX,  M,  — l.,l.) 

CALL  DETER ( CAX, M,  CNX ) 

CALL  DETER ( CBX, M,  CSX) 

RATX=REAL ( CNX /CSX ) 

WRITE(6,  3)  RATX, THRX 

FORMAT (/,  IX,  "TENS: RATIO*" , F9. 4, IX, "THRESHOLD*",  F7. 4) 
TWENTIES  DETECT 
IF( J. EQ. 1 )  GOTO  4 

CALL  PRODV ( CVXX, M,  CVXX,  M, CBXX, WXX ) 

CALL  CLRV(CVXX,  M) 

CALL  ADDM (CBXX, CAX X,M, -1., 1. ) 

CALL  DETER(CAXX, M,  CNXX) 

CALL  DETER ( CBXX,  M,  CSX X) 

R AT  XX=REAL (CNXX /CSX  X ) 

CALL  CLRM ( CAX  X , M, M ) 

CALL  CLRM ( CBX  X , M, M ) 

WRITE(6,5)  RATXX, THRXX 

FORMAT ( /,  IX,  "TWENTIES: RATIO*", F9. 4, IX, "THRESHOLD*",  F7. 4) 
FIFTIES  DETECT 


IF ( I*J. EQ. 3. OR. J*I . EQ. 10)  GOTO  8 
GOTO  1 

CALL  PRODV  ( CVL,  M,  CVL,  M,  CBL,  WL ) 

CALL  CLRV( CVL, M) 

CALL  ADDM (CBL, CAL,  M,  -1 . ,  1. ) 

CALL  DETER ( CAL, M,  CNL) 

CALL  DETER ( CBL,  M,  CSL ) 

RATL=REAL ( CNL/CSL ) 

CALL  CLRM ( CAL, M,  M) 

CALL  CLRM ( CBL, M,  M) 

WRITE (6, 6)  RATL, THRL 

FORMAT ( /, IX,  "FIFTIES: RATIO-', F9. 4,  IX,  "THRESHOLD*",  F7.  4 > 
CONTINUE 

HUNDREDS  DETECT 


f rnnti nimH  rm  npyt  nanp^ 


CALL  PRODV  ( CVC,  M,  CVC,  M,  CBC,  WC ) 

CALL  CLRV(CVC,M) 

CALL  ADDM ( CBC, CAC,  M,  - 1 . ,  1 . ) 

CALL  DETER ( CAC, M, CNC ) 

CALL  DETER ( CBC, M,  CSC ) 

RATC=REAL ( CNC/CSC ) 

CALL  CLRM ( CAC, M, M) 

CALL  CLRM ( CBC, M, M) 

WRITE (6,  7)  RATC,  THRC 

FORMAT (/,  IX,  "HUNDREDS: RATIO*", F9. 4,  IX, "THRESHOLD*",  F7. 4 > 

STOP 

END 


F/6  12/ 1 


A0-AQ91  »54 


unclassified 
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■ 


o  o  o 


I 


PROGRAM  P-6  (continued) 


SUBROUTINES 


» 


I 


) 


> 


I 


SUBROUT I NE  DETER ( CMAT, M, CDET ) 

IMPLICIT  COMPLEX <C) 

DIMENSION  CMAT ( M,  M  ) 

IMAX=M-1 
CDET=< l.,0. ) 

DO  1  1=1,  IMAX 
CPIV=CMAT ( I , I ) 

CDET=CDET*CPIV 
JMIN=I+1 
DO  1  J=JMIN, M 
DO  1  J1=JMIN,  M 

1  CMAT ( J,  J1 ) =CMAT  <  J,  J1 ) -CMAT  < I , J1 ) »CMAT ( J,  I ) /CPI V 

CDET=CDET«CMAT  <  M, M  > 

RETURN 

END 

C 

SUBROUT I NE  PRODV  <  CV 1 ,  N1 , CVS, N2,  CMAT,  W ) 

IMPLICIT  COMPLEX <C) 

DIMENSION  CV1 < N1 ) , CV2( N2 > , CMAT ( N1 ,  N2> 

DO  1  1  =  1,  N1 
DO  1  J= 1 ,  N2 

1  CMAT ( I,  J) =W*CV 1(1) KCONJG ( CV2 ( J ) ) 

RETURN 

END 

C 

SUBROUT I NE  ADDM ( CMAT 1 , CM ATS,  M,  W 1 ,  W2 > 

IMPLICIT  COMPLEX (C) 

DIMENSION  CMAT 1 <M,  M ) ,  CMAT2 ( M,  M) 

DO  1  1  =  1,  M 
'  DO  1  J= 1 ,  M 

CMAT  1  ( I ,  J )  =W  1  **CMAT  1  ( I ,  J )  +W2*CM ATS  ( I ,  J ) 

1  CONTINUE 

RETURN 
END 
C 

SUBROUT I NE  CLRM ( CMAT, M,  N  > 

IMPLICIT  COMPLEX (C) 

DIMENSION  CMAT (M, N) 

DO  1  1=1, M 
DO  1  J=l,  N 

1  CMAT (I,J)  =  (0.,0.  ) 

RETURN 

END 

C 

SUBROUT I NE  CLRV ( CV,  M ) 

IMPLICIT  COMPLEX <C) 

DIMENSION  CV(M) 

DO  1  1=1, M 
1  CV(I)=(0.,0. ) 

RETURN 

END 


I 


(continued  on  next  page) 


PROGRAM  P-6  (continued) 

C 

SUBROUTINE  ADDV <CV1,  CV2,  M,  W1 r  W2 ) 
IMPLICIT  COMPLEX (C) 

DIMENSION  CV1 (M),  CV2<M) 

DO  1  1=1, M 

1  CV1<I>=W1*CV1(I)+W2*CV2(I> 

RETURN 

END 

/* 

//LKED.SYSLIB  DD 

/ /  DD  DSN=C59 1 . SUBL I B,  D I SP=SHR 

/* 


OUTPUT  FROM  PROGRAM  P-6  AS  LISTED,  SHOWING  DETECTIONS  (Pc-=  . 
(5  sensors,  S  =  .3  or  SNR  =  -13.5  dB)  PA 


TENS: RATIO*  3.0335  THRESHOLD*  5.S4S6 

TENS: RATIO*  1.6548  THRESHOLD*  5.8486 

TWENTIES: RATIO*  1.6419  THRESHOLD*  1.9930 
TENS: RATIO*  3.5713  THRESHOLD*  5.8486 

TENS: RATIO*  3.8445  THRESHOLD*  5.8486 

TWENTIES: RATIO*  1.9891  THRESHOLD*  1.9930 

TENS: RATIO*  3.0560  THRESHOLD*  5.8486 

FIFTIES: RATIO*  1.5326  THRESHOLD*  1.2805 

TENS: RATIO*  3.0316  THRESHOLD*  5.8486 

—>  T WENT I ES : RAT I 0=  2.1577  THRESHOLD*  1 . 9930 

TENS: RATIO*  1.5716  THRESHOLD*  5.8486 

TENSIRATIO*  3.7013  THRESHOLD*  5.8486 

TWENTIES: RAT 10=  1.3103  THRESHOLD*  1.9930 

TENS: RATIO*  1.7447  THRESHOLD*  5.8486 

TENS: RAT 10=  1.7189  THRESHOLD*  5.8486 

TWENTIES: RATIO*  1.2670  THRESHOLD*  1.9930 

1.1324  THRESHOLD*  1.2805 
1.2309  THRESHOLD*  1.1271 


FIFTIES: RATIO* 
— *  HUNDREDS: RATIO* 


APPENDIX  A 


DISTRIBUTION  OF  MULTIVARIATE  TEST  STATISTIC  FOR 
COMPLEX  DATA 


(A)  Following  the  approach  of  Anderson,  we  note  that  the  sample 
mean  vector  and  sample  covariance  are  independent.  The  original  data 
vector  samples  are 

Xfc  =  +  JVfc.  (m  x  1).  ( A- 1 ) 

v/here  samples  are  assumed  independent  and  identically  distributed  as 
multivariate  normal,  denoted 

~  N(a,z),  v^  ~N(b,2),  independent).  (A-2) 

A 

The  sample  covariance  £  is  such  that 
2nZ  =  XX*  -  X0X0*  =  A 

n  n 

=  "  n-^*  ’  £  =  n  X)  *A'3^ 

k=l  k=l 


The  data  matrix  can  also  be  represented  in  a  different  coordinate  system; 
let  B  be  an  orthogonal,  (n*n)  matrix: 


X  =  ZB 1 ,  BB'  =  In;  B=(br  b2 . b^.  (real)  (A-4) 

That  is,  there  is  another  set  of  vectors  such  that 

n 

Ik  ’  %  m  X)  biiA 

i=l 


and  these  vectors  are  also  Independent.  Let  the  nth  column  of  B  have 
all  components  equal  to  1/v'n  .  Then 

N(/n  i) 


t/n  G 


~  N(0,  s),  Mn 


(A— 5) 


A  =  ZB'BZ*  -  n  v  v*  »  Z  Z*  -  z_  z  * 


-n  -n 


4  4*  -  V  ■  Z/ 


2  2  * 
4  4  • 


(A-6) 


Therefore  ju  is  independent  of  A  because  z  is  independent  of  the  other 


(B)  Using  these  new  vectors,  the  test  statistic  is 

h  v)'1  h 

4' 

Now  the  {z^}  are  transformed  by  diagonalization: 

^  where  DED'  =  Im 

so  that  n  £*  A"1  £  =  tr  j( ^  ^  j^*)’1  ^ 


(A-7) 


(A-8) 

(A-9) 


where  the  have  identity  covariance  matrices. 

(C)  A  third  transformation  is  performed  using  a  unitary  matrix 
(complex  orthogonal)  Q  such  that  QQ*  =  I  and  resulting  in  the  data 


vectors 


*k  = 


with  the  constraint 


(A-10) 


=  0^  =  (|*J,  0,  0 . 0)'.  (A-ll) 

This  makes  the  matrix  ^  in  (A-9)  become  a  matrix  with  only  one 

nonzero  component,  causing  the  trace  to  be 


n£*  A  *  £  = 


fa  V1 

2-r  *k  V) 

tfc=l  /  _ 


The  element  b1*  of  the  inverse  of  B  *  ^  w^ 


is  given  by 


jTT  =  bll  ^1)  B22  ^(1) 


where  B  has  been  partitioned  into 


”bll 

_bd 


L(D  22  J  * 

Therefore  the  test  statistic  is  given  by 


A  'tr  n  •  _ 

•njj  A  a  - 


bll  "  ^fl)B22  ^{1) 


(A-12) 


(A-13) 


(A- 14) 


(A-15) 


The  numerator  can  be  seen  as  a  noncentral  chi -squared  variate 
IwJ2  ~  X2(2m;X) 
with  noncentrality  parameter 

*  E{*n*}  EK}  ’  EM  D'0* .«  D  EW 

=  n  £*  D*  D  £  =  n  u*  x~*  u 


(A- 16) 


By  Theorem  4.3.3  of  Anderson  the  denominator  of  (A-15)  is  seen  to  be 
(central)  chi -squared,  independent  of  the  numerator: 

l/bH  ~  xZ(2n-2m).  (A-17) 

Thus  the  ratio  is  an  F-statistic,  denoted 


!2m;x)  _  m 


x*(2n-2m) 


^9m  9 


n-m  2m,2(n-m)' 


•  *  A-l  - 

n±L  A  £ 


(A-18) 


r 


Appendix  B 

Moments  of  the  Random  Signal  Test  Statistic 


(A  completely  analogous  derivation  for  real  data  is  given  by  Anderson 

[2].) 

For  n  samples  from  an  m-dimensional  distribution  of  zero-mean,  complex, 
jointly  Gaussian  variables  with  covariance  matrix  z  the  sample  covariance 


A  =  2n  z  - 

k=l 

has  the  Wishart  distribution  [2,  3] 


P(A;n)  *  yz.n)  |A|n'm  etr 


7  e"1a 


where 


The  moments  therefore  of 

2  *  |A0I/|A|  ,  A0  I  d1ag(an, 
are  given  by 

E{zh}=  /dA|A0|h  |ft|'h  p(ftsn) 


4  r1* 


=/dA  Km(E,n)|A0|1’  |A|n_m_h  etr 

Kjz,  n)  r  h 
=  yz,  n-h)  J  dA  lAol  p(A »n-h) . 

Integrating  first  over  the  off-diagonal  elements  yields 
Elzf1*  Kjz,  "ri-h)  f  dA0  lA0»h  p0(A0in_h)* 


(B-l) 


yz.n)'1:  2mn  Trn,(m-1)/2  r(n)r(n-l). ..  (n-m+l)|z|n.  (B-2) 


(B-3) 


(B-4) 


(B-5) 


where  pQ(  )  is  the  joint  pdf  of  the  diagonal  elements  of  A.  For  z 
diagonal  this  density  is  extremely  complex  in  expression  [9,  10]. 
for 


z 


°22'  *  ‘ 


mm 


)• 


these  elements  are  independently  distributed  as 


siin'1'h  exp{-an/2onl 


and  the  integration  is  simple,  producing 


„  *Sn(z0*n)  _  i 

i=l 


Kihi-n-h) 

Kl(51VT) 


_  r(n-h)r(n-h-l) . . .r(n-h-m+l) 
r(n)r(n-l)...r(n-m+l) 


r(n) 

r(n-h) 


This  expression  may  be 


further  developed  as 


r(n-h-i+l)r(n) 

r(n-i+l)r(n-h) 


B(n-h-l+l.  1-1)  f 
B(n-i+l,  i-1) 


where  the  beta  function  is 


non- 

However , 

(B-6) 

(B-7) 


CM) 

(B-9) 


B(a,  b)  =  liiilM. 

r(a+b) 


(B-10) 


Appendix  C 

Effect  of  Correlation  and  Phase 
on  Multivariate  Detection  Parameter 


f: 


From  (4-10),  the  multivariate  detection  parameter  (generalized 
SNR)  is  given  by 


m  m 


AM  =  n 


air  s^Sr  cos(e..  -  er). 


(C-l) 


i=l  val 


where  jo^rJ  *  is  the  inverse  of  the  intersensor  noise  covariance 
matrix  and  the  |s.j ,  e^j  are  the  amplitudes  and  phases  of  the  signals 
arriving  at  the  sensors  (or  magnitude  and  phase  of  the  DFT  bin  output). 

If  the  phases  are  independent  and  uniformly  distributed  on  some 
interval  (-A,  a),  then 

E{cos(ei  -  er)J  =  .  (C-2) 


Thus  we  have 


EeM 


"Z  SiC"+  "  Z  SiSr  °ir* 

i  A  “ 


(C-3) 


and  the  second  term  vanishes  for  A  =  kw.  A  random  sensor  placement  therefore, 
such  as  shown  in  Figure  C-l,  could  cause  the  second  term  to  drop  out  regard¬ 
less  of  the  covariance  matrix  structure. 

If  the  covariance  matrix  E  can  be  represented  by  the  case  in  which 
the  intersensor  correlations  are  all  equal,  then 

£  *  “  DoR  °o*  (c-4) 


using 


Do  =  diag(o^. 


(C-5) 


and 


Sr}  • 


*i1  =  *•  1  *  r 

*ir  -  p,  i  J1  r 


(C-6) 


For  this  special  case 


E_1  =  D_1  R-1  D"1 


(C-7) 


and 


ir 


RaZ2  ,  i  =  r 
a/aiar,  i  /  r 


(C-8) 


where 


B  = 


l+(m-2)p 
(l-p)[l+(m-l)p] 


and 


(C-9a) 


°  =  (r-pTCl^m-Dp]  • 

In  this  instance  the  average  detection  parameter  has  the  form 

MV”}  =  nB  *  sl/ol  *  ne(STi)*(*  V8«)2 


(C-9b) 


2nm  SNR 


avg 


[0  ♦  (m-l)a(^)2]. 


(C-10) 


Appendix  D 

Distribution  of  Sample  Correlation  Coefficient  (Squared) 


The  elements  of  the-sample  covariance  matrix  2nt  =  A  for  m=2  have 
the  Wishart  distribution 

i  n-2 


Pj (A)  =  K(Z,n)  |A|n_<:  etr 


where 


K_1(z,n)  =  2cn  ir | z | nr(n)r(n-l) . 


=  o2n  _ i ^  t  n 


(D-l) 

(D-2) 


Explicitly  (D-l)  is  written 


n-2  | 

1  1 

exp 

"2,2 

'  1 

1-p  l 

Pl(A)  ’  4n,[o5o|(l-p2)]n'  r(n)r(n-l) 

Let  a  transformation  of  variables  be  defined  by 


*11  9  a12R  .  a22, 

~r  *  2po'^r +  ~r\ 

Oj  1  £.02 


(D-3) 


2  la12l  /aHa22 


ucos0  =  a 
usinQ  =  a 


=  a 


12R 

121 

22 


|J|=  JS- 


vz 


The  joint  pdf  of  the  new  variables  is 

2(n-2)+3,.  ,n-2 

P2(z,u,e,v)  =  K  ^  ■  ■  ■  -  -  exp 


(D-4) 


ML  +  ,__p_u_cos9  _  , 

vz  (\ 2T  1 

al°2  ^”P  ' 


In  the  next  steps  the  unwanted  variables  are  integrated  out: 
2ir 


(D-5) 


/. 


,  .  2(n-2)+3  /,  %n-2 

P2(z,u,e,v)do  =  2irK  — 

0  2 


pu 


°1°2(1”P  ^ 


BU 

vz 


-aV 


=  p3(z,u,v) 
(D-6) 


00 


fl±A"'2e-«'>f  dx  x«-l  i0[r^r]  ,-««/« 
z  J  0 


=  r(n)irK 


-aV 


=  P4(z»  v) 
(D-7) 


oo 


TrK(l-z)n~2  r(n)r(n) 


l(n,n;1;  4fSa) 


-  rxMy  (1-p2)"  (1-z)n‘2  2Fl(n’n;1;p2z)  =  p5(z) 

(D-8) 


The  complementary  cumulative  distribution  function  for  z  is 
computed  as  follows: 


2k 

klkT  B^k+1*  n~l) 1  i_n (n-l  »  k+D 

(D-9) 


or 


Qz(n)  ■  (l-»2)" 


2k 

V^k 


Ij_r)(n-1 ,  k+1) 


(D-10) 
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Glossary  of  Notation 


a^ 

A 

b 

Br(x) 

B(a,b) 

det(B) 

OFT 

E{x} 

etr(B) 

F 

F 

VI, v2 

jFj(a;b;x) 

2^  j(a  **5»c 

H0» 

I 

Ix(a,b) 

j 

A(X) 

X 

m 

MDS 

ML 

MPG 

1L 

n 

pdf 


real  part  of  mean  vector 
sample  covariance  matrix  (2nr) 
imaginary  part  of  mean  vector 
Bernoulli  polynomial 
beta  function 

determinant  of  the  matrix  B 
discrete  Fourier  transform 
mathematical  expectation  of  quantity  x 
exp{tr(B)} 

matrix  of  OFT  components 

F-statistic  with  degrees  of  freedom  vj  and  v2 

confluent  hypergeoroe trie  function 

Gaussian  hypergeometric  function 

hypotheses 

identity  matrix 

Pearson's  incomplete  beta  function 
/IT 

likelihood  ratio  value 

noncentrality  parameter,  detection  parameter 

number  of  sensors 

minimum  detectable  signal 

maximum  likelihood 

multivariate  processing  gain 

mean  vector 

number  of  samples 

probability  density  function 

probability  of  detection 


probability  of  false  alarm 

probability  of  event  A 

probability  integral 

receiver  operating  characteristics 

correlation  coefficient 

covariance  matrix 

signal -to-noise  ratio 

trace  of  matrix  B 

real  part  of  data  vector,  matrix 

imaginary  part  of  data  vector,  matrix 

data  vector,  matrix 

statistic,  test  value 

estimated  value 

column  vector 

transpose  of  vector,  matrix 
conjugate  transpose  of  vector,  matrix 
determinant  of  matrix  B 
magnitude  of  vector 


